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Preface 


Nuclear thermal hydraulics is the application of thermofluid mechanics within the nuclear indus- 
try. Thermal hydraulic analysis is an important tool in addressing the global challenge to reduce 
the cost of advanced nuclear technologies. An improved predictive capability and understanding 
supports the development, optimisation and safety substantiation of nuclear power plants. 


This document is part of Nuclear Heat Transfer and Passive Cooling: Technical Volumes and Case 
Studies, a set of six technical volumes and four case studies providing information and guidance 
on aspects of nuclear thermal hydraulic analysis. This document set has been delivered by Frazer- 
Nash Consultancy, with support from a number of academic and industrial partners, as part of 
the UK Government Nuclear Innovation Programme: Advanced Reactor Design, funded by the 
Department for Business, Energy and Industrial Strategy (BEIS). 


Each technical volume outlines the technical challenges, latest analysis methods and future direc- 
tion for a specific area of nuclear thermal hydraulics. The case studies illustrate the use of a subset 
of these methods in representative nuclear industry examples. The document set is designed for 
technical users with some prior knowledge of thermofluid mechanics, who wish to know more about 
nuclear thermal hydraulics. 


The work promotes a consistent methodology for thermal hydraulic analysis of single-phase heat 
transfer and passive cooling, to inform the link between academic research and end-user needs, 
and to provide a high-quality, peer-reviewed document set suitable for use across the nuclear 
industry. 


The document set is not intended to be exhaustive or provide a set of standard engineering ‘guide- 
lines’ and it is strongly recommended that nuclear thermal hydraulic analyses are undertaken by 
Suitably Qualified and Experienced Personnel. 


The first edition of this document set has been authored by Frazer-Nash Consultancy, with the 
support of the individuals and organisations noted in each. Please acknowledge these documents 
in any work where they are used: 


Frazer-Nash Consultancy (2021) Nuclear Heat Transfer and Passive Cooling, 
Volume 2: Convection, Radiation and Conjugate Heat Transfer. 
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1.1 


Introduction 


Conjugate Heat Transfer (CHT) is the coupling of heat conduction in solids with convective and 
radiative heat transfer in neighbouring fluids. The ability to predict the variation in temperature of 
solid components, such as pipework, vessels and structures, as a result of heating or cooling by a 
fluid is important for all Nuclear Power Plants (NPPs). 


The objective of this technical volume is to describe how to model convection, thermal radiation 
and CHT, which are fundamental to the transfer of heat out of a reactor particularly under passive 
cooling scenarios. It is part of a set of documents that are intended to describe good practice and 
summarise the current state-of-the-art, and it is recommended that Volume 1 (Introduction to the 
Technical Volumes and Case Studies) is reviewed for context. 


The benefit of CHT analysis and its impact on reactor design is discussed in this section before 
providing a detailed review of the modes of heat transfer and current analysis methods. 


Conjugate Heat Transfer 


CHT analysis is required when it is not sufficient to simply consider heat transport in the fluid, with 
an idealised thermal boundary condition applied at the fluid and wall interface. In these cases, the 
fully coupled problem, including heat conduction in the solid wall adjacent to the fluid, must be 
considered. Situations where this is necessary include: 


¢ When the fluid temperature varies during a transient, the conduction behaviour and thermal 
mass' (stored energy) inside the vessel walls may have a significant impact on the heat 
transfer to the fluid and surroundings, and therefore the system performance as a whole. 


* Conduction and the thermal diffusivity of the solids may have a significant impact on the flow 
behaviour if the temperature field inside a fluid and the surrounding vessel or pipework is 
complex, unsteady or three-dimensional. 


+ If the fluid temperature varies adjacent to a wall (for example due to unsteady poorly mixed 
flow or a moving free surface) the resulting temperature variations within the solid may cause 
fluctuating thermal stresses that can lead to thermal fatigue. 


Therefore, in many circumstances it is necessary to properly account for the temperature variations 
within solids and fluids flowing around or through solid components, and the thermal coupling 
between the two. This can be achieved using a CHT approach, which predicts the temperatures 
within a fluid and the adjacent solids within the same analysis or coupled together. 

Thermal mass (or heat capacity) refers to the product of mass and specific heat capacity of a body (mcp) and is a measure 
of a body’s ability to absorb energy. This is sometimes referred to as ‘thermal inertia’, but this is a separate property 


(,/kpcp) used in geotechnical or planetary science but rarely in mechanical engineering. Thermal diffusivity (a = k/pcp) 
is a more relevant related quantity, which has significance in the transient conduction equation. 
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Introduction 


Impact on Reactor Design 


The accurate prediction of fluid and solid temperatures is used to support performance assess- 
ments and structural integrity calculations as part of the design and safety of all NPPs. 


Performance Assessments 


Detailed predictions of the fluid and solid temperatures within the primary circuit are required under 
normal operation, start-up and shut-down transients, as well as during fault scenarios. These are 
used to determine the operating envelope and safety margins for the reactor design. This is partic- 
ularly relevant to passive cooling systems, where the development of thermal gradients is essential 
to generate the buoyancy forces needed to drive the flow through the system. For example, CHT 
analysis is used to: 


« Ensure that the maximum temperatures of key components remain within allowable limits 
(e.g. fuel cladding, fuel pellets, reactor pressure vessel and primary circuit components). 


* Determine the rate of cooldown of the reactor under normal shut-down and fault scenarios 
and demonstrate the effectiveness of passive cooling systems, as a prolonged cooling period 
is required due to the large thermal mass and ongoing nuclear decay reactions. 


« Predict the temperatures or heat flux where the conduction paths are complex/subtle, such 
as some decay heat removal systems or cross-core conduction under passive cooling. 


Structural Integrity Calculations 


CHT analysis is a key input into structural integrity calculations (Section 2.4) in order to demon- 
strate the integrity of Structures, Systems and Components (SSCs) and predict component life. 
Failures in components may be caused by situations such as large thermal shocks causing fracture 
(e.g. emergency cooling scenarios), temperature fluctuations over a long period causing fatigue or 
temperatures that lead to enhanced corrosion. Some examples relevant to NPP design include: 


* Thermo-mechanical studies are required to calculate component temperatures, as internal 
stored energy, heat generation and conductivity induce stresses in components, welds and 
joints, especially if their state changes rapidly. 


* Temperature fluctuations induced in solid materials through turbulent mixing processes may 
lead to thermal fatigue (cyclic thermal stresses and eventually to fatigue and cracking). For 
example, in 1998 a through-wall crack was discovered downstream of a T-junction (Kuhn 
et al., 2010) at the Civaux-1 plant in France only five months after finishing construction, 
resulting in pipe rupture and release of radioactive steam into the reactor building. Therefore, 
predicting the impact, location, amplitude and frequency of those temperature fluctuations is 
crucial for the lifetime management of components used in NPPs. 


« Excessive thermal stress caused by temperature differences in the fluid (e.g. thermal strip- 
ing, stratification) can lead to thermal fatigue damage and cracks in feedwater nozzles, high 
pressure Safety injection lines and residual heat removal lines (Jo et al., 2003). For example, 
Pressurised Thermal Shock (PTS) is a key reactor safety issue that may impact the lifespan 
(and subsequent economics) of the plant due to stresses in the Reactor Pressure Vessel 
(RPV) downcomer caused by cold water injection (IAEA, 2010). 
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2 Technical Context 


This section considers some of the key technical aspects of heat transfer with a focus on nuclear 
applications. The modes of heat transfer are introduced, followed by a discussion on surface mod- 
ifications, material properties, thermal fatigue and key modelling challenges. 


2.1 The Modes of Heat Transfer 


Heat transfer is the movement of thermal energy due to a temperature difference. There are three 
modes of heat transfer, which are illustrated in Figure 2.1: 


Conduction: In conduction, energy is transferred from more energetic particles to less energetic 
particles by interactions between them (e.g. direct interaction between atoms or molecules, 
movement of free electrons in metals). Conduction is generally the main mode of heat trans- 
fer within compact (non-porous) solids. 


Convection: In convection, energy is transferred by the combination of bulk motion of a fluid (ad- 
vection) and by diffusion (resulting from random molecular or turbulent motion). Convection 
is generally the most significant mode of heat transfer between solids and fluids at moderate 
temperature differences. 


Thermal radiation: In thermal radiation, energy is transferred by the electromagnetic radiation 
emitted from matter, so no medium is required and heat transfer can occur in a vacuum. 
Solids, liquids and gases can all participate in thermal radiation heat transfer. Thermal radi- 
ation is often most significant at high temperatures. 


Warm Solid Insulation Coal 
Working Fluid Surroundings 
—» Conduction 
wo _— 
—a4 — —> Convection 
—_ __ _> _ —~y» Thermal radiation 


Steady-State 
Temperature Trend 


Figure 2.1: Example of three modes of heat transfer at an insulated wall. 


One, two or all three of these modes may be important in a given situation. These modes are 
discussed further in the following sections. More general detail and explanation on heat transfer is 
available in Rogers and Mayhew (1992), Incropera et al. (2011) or Rohsenow et al. (1998). 


1 In addition, Lienhard and Lienhard (2020) is free to download from web.mit.edu/lienhard/www/ahtt.html. 
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Conduction 


The rate of heat transfer per unit area or heat flux (q) through a substance due to conduction is 
proportional to the gradient in temperature (7) according to Fourier’s law: 


dT(x) 


q=-—kVT_ orinonedimension, qx = —k, 7d 
x 


Where k is the thermal conductivity (assumed constant in this equation), a property of the material 
(Section 2.3). Although this is usually expressed as a scalar, it may be necessary to use a ten- 
sorial description for anisotropic materials, such as composites (Section 2.3) or ‘effective’ material 
properties in porous regions (Section 3.4.3.3). 


Thermal Resistance 


In heat transfer analysis it is often useful to separate the details of the geometry and material 
properties (which are often fixed in a given situation) from the temperatures (which may vary). This 
can be achieved for steady-state analysis using the concept of thermal resistances. In the simple 
plane wall in Figure 2.1, k is assumed to be a constant and the temperature gradient through the 
wall is therefore linear. Applying Fourier’s law to such a wall, the rate of heat transfer through a 
given frontal area (A, normal to the heat flux) for wall thickness (L) is: 
Q= kan — zaAT where Rip = - 

Where AT is the temperature difference between the two surfaces of the wall and R,, is the thermal 
resistance due to conduction through the wall (which is large for a well insulating wall). Different 
basic geometries have different formulae for this thermal resistance (Section 2.1.2.4 presents an 
example for radial conduction in a cylindrical tube or pipe). 


Thermal resistances are convenient because they can be used in a broadly similar manner to 
electrical resistances (i.e. by adding resistances in series and reciprocals of resistances in parallel). 
Section 2.1.2 and Section 2.1.3 develop thermal resistances further and include convection and 
thermal radiation. All of the thermal resistances occurring at the wall (for the three modes of heat 
transfer in Figure 2.1) can therefore be combined to form an overall thermal resistance Rp, oa- 


This overall thermal resistance can be used to calculate the widely used ‘overall heat transfer 
coefficient’, ‘overall thermal transmittance’, ‘U-value’ or U (with units Wm? K~!). However, it is 
important to note that approaches to defining R:, vary. In this volume R:; includes A (its units are 
K W-*) because for geometries where the frontal area changes along the direction of the heat 
flux (like radial conduction in a pipe) this is clearer. Therefore, VA = 1/R¢h,oa, where A is some 
representative area, the basis of which should be clearly stated*, see Section 2.1.2.4. 


In some references (particularly associated with the construction industry, which usually works with plane walls where the 
value of A is clear), Rin may not include A (i.e. have units K m? W-!) and accordingly U would equal 1/R¢n. 
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2.1.1.2 


Technical Context 


Contact Thermal Resistance: Two solids that are touching are rarely in perfect thermal contact. 
Instead, the roughness of their surfaces and possibly small distortions in their geometry can cause 
narrow voids at the interface between the solids. These voids are likely to cause an additional 
thermal resistance. It may be appropriate to assess this contact thermal resistance by considering 
two resistances in parallel: 


1. Conduction through the (probably small) patches where there is good thermal contact. 


2. Conduction (and possibly also convection and thermal radiation) through the voids at the in- 
terface (which are likely to be filled with air or the surrounding fluid) based on a characteristic 
gap width. 


The contact resistance may also be affected by forces acting over the interface, wear occurring 
in service, the cleanliness of the surfaces, and any other materials present at the interface. For 
more detailed assessments, experimentally measured contact resistances are available in wider 
literature (Rohsenow ef al., 1998). 


Unsteady Conduction 


Conduction through solids during transients is a key aspect of CHT, and this heat transfer is a 
function of temperature gradients in the solid that are likely to change over time. It is therefore 
useful to develop an understanding of the behaviour of these temperature gradients where a solid 
is in contact with a fluid, for which a key non-dimensional group is: 


Biot number (Bi = hL/k,), a ratio of convective and conductive heat transfer effectiveness 


Where h is the convective Heat Transfer Coefficient (HTC) (Section 2.1.2), L is a characteristic 
length scale (e.g. half the thickness of a plate or the volume to surface area ratio for a more 
general geometry) and k, is the solid thermal conductivity. Estimating Bi using a HTC correlation 
is often useful in the early stages of a CHT analysis: 


¢ If Bi < 0.1, conduction within the solid is much more effective than convection to the fluid, so 
temperature gradients within the solid are likely to remain small. It may therefore be appropri- 
ate to ignore spatial temperature variations within the solid and use a single mass-average 
‘bulk’ or ‘mean’ temperature to predict heat transfer through time. This is sometimes referred 
to as a ‘lumped capacitance’ approach. 


« If the above criterion is not met, the surface of the solid is likely to respond to the fluid tem- 
perature more rapidly than the interior of the solid. As Bi increases, larger thermal gradients 
would be expected within the solid, so both spatial and temporal variations in temperature 
within the solid may be important for predicting heat transfer. 


These two scenarios above are illustrated in Figure 2.2. 


Another non-dimensional group that is often encountered in CHT analysis is the Fourier number 
(Fo = a,t/L?), a non-dimensional timescale for conduction analysis (Rogers and Mayhew, 1992 
and Incropera et a/., 2011). In the limiting case where the Biot number approaches zero, a ‘lumped 
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Temperature 
Trend (later) 


Surrounding Fluid 
(Warming) 


Temperature 
Trend (initial) 
Small Bi Large Bi 


Figure 2.2: Illustration of temperature gradients in a solid surrounded by a warming fluid 


at two points in time, for small and large values of B/. 


capacitance’ approach can be used to calculate the rate of heating or cooling of the solid: 


Bie = e-BiFo _ -(we)t : eee AS: In @ 


6;  T;—Too a hA, 0; 


Where 7; is the initial temperature of a solid surrounded by a fluid at temperature 7,,., and T is the 
temperature of the solid at time t. 


Convection 


Fluid motion is generally a key aspect of convective heat transfer, and a distinction is made between 
external and internal flows. Engineering judgement may be required to understand whether a flow 
is best considered as internal or external flow. For example, for a small cylinder crossing a large 
duct containing essentially uniform flow, the flow around the cylinder may closely resemble external 
flow, even though the flow down the duct as a whole is an internal flow. 


External Flow: In this case, the flow as a whole is considered to be unconstrained (e.g. the flow 
over an aeroplane wing). There is therefore normally a clear ‘free-stream flow’ that is unaffected 
by surfaces, and a distinct boundary layer flow commonly develops adjacent to surfaces. The local 
heat flux per unit area (q) at a surface due to convection can be considered to be proportional to 
the difference between the temperature of the surface (7,,) and some representative free-stream 
fluid static temperature (T;,.0) for low speed flows?, according to Newton’s law of cooling: 


q = h(Tw — Ts,co) 


Where h is the convective HTC, which depends on the details of the flow near the surface (Sec- 
tion 2.1.2.1). 


Internal Flow: In this case, the flow as a whole is constrained by surfaces (e.g. the flow in a 
pipe). Significant gradients may therefore occur throughout the flow field, and there is often no 
clear free-stream flow. There is no clear 7, .. that can be used for Newton’s law of cooling, so it is 
therefore necessary to define a ‘bulk’ or ‘mean’ static temperature (7, or T,) that is characteristic 
of the flow and can be used instead of T; ... This bulk temperature is commonly based on a quasi- 


For high speed (high Ma) flows, the adiabatic wall or recovery temperature T;,, is typically used instead of Ts,o0, to account 
for kinetic energy in the flow being dissipated as heat in the boundary layer. 
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mass-average across the flow’: 


_ oe T;dW 1 


s—= "Sf aay = tars UTsdAcs 
lig dw WwW Acs ‘ 


The local heat flux per unit area (q) at a surface due to convection in internal flow can then be 
calculated using a HTC in a similar manner to low speed external flows: 


q= h(Tw = Ts) 


It is noted that while 7, is generally constant in low speed external flows, in internal flows T. 
varies as the flow is heated or cooled (e.g. along the flow in a pipe). 


Near Wall Flows and Convection 


Immediately next to a surface (or ‘at the wall’), the velocity and temperature of the flow is normally 
considered to be the same as those at the wall itself, so heat transfer at the fluid-solid interface 
is entirely due to heat conduction. Further into the fluid, the wall-parallel fluid motion enhances 
the transfer of thermal energy, which makes the temperature profile in the fluid layers away from 
the wall flatter and the temperature gradient at the wall stronger. Comparing Fourier’s law (Sec- 
tion 2.1.1) and Newton’s law of cooling (above) indicates that the HTC depends on this temperature 
gradient at the wall. This section introduces some aspects of flows that can alter these gradients, 
and hence affect heat transfer. 


As noted above, external flows often feature a free-stream flow with momentum and thermal bound- 
ary layers between this flow and the wall. Across these two boundary layers, velocity and tempera- 
ture respectively, start to change from their wall values to asymptotically approach their free-stream 
values. In the case of the velocity boundary layer, the rate of growth depends on the balance of 
the rate at which the wall-parallel momentum is convected along the wall-parallel direction and the 
rate at which it is transported in the wall-normal direction through molecular and turbulent mixing. 
The rate of growth of the thermal boundary layer is determined by the corresponding balance in 
the transport of thermal energy (Figure 2.3). 


The temperature and velocity gradients (and hence HTC) are well understood for simple boundary 
layers, such as flat plates or fully developed pipe flow, but are strongly affected by factors such as: 


¢ The nature of the free-stream flow. Aspects such as local acceleration and deceleration, 
turbulence levels, vorticity and flow impingement on the surface. 


« Whether the flow is laminar, turbulent or transitional. Turbulence enhances the mixing of 
momentum and heat (across velocity and temperature gradients respectively). 

¢ The behaviour of the boundary layer itself. Aside from turbulence, features like separation 
have a significant impact. 

« The properties of the fluid. These have a significant impact on the transport of momentum 


and heat. 


Strictly, only extensive properties (that define a property per unit mass) can be mass-averaged. In this case, the relevant 
extensive property is enthalpy, so strictly this equation is only valid for a perfect gas (Volume 3, Section 2.3.1) or incom- 
pressible liquid with a constant specific heat, since in these cases a change in enthalpy equals cp AT; and cp can be 
cancelled from the equation. Further details on averaging non-uniform flows is available in Cumpsty and Horlock (2005). 
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Momentum 
= S S 
= aes Hey 
L EY Ey, 
Freestream Momentum Thermal |Z LL Ez 
Thermal 
ie > > 
i jo 
ey, ey. 
of Ez Bz 
(a) External Flow (b) Internal Pipe Flow 


Figure 2.3: Illustration of momentum and thermal boundary layers. 


* The detailed geometry of the surface. Roughness and surface features (Section 2.2) can 
have a significant impact on the local flow field near the surface. 


The behaviour of boundary layers is a significant topic in its own right, and a detailed discussion 
is provided by Schlichting and Gersten (2017). Turbulence and transition both affect heat transfer 
and are particularly significant issues for passive cooling systems; further details are provided in 
Volume 3 (Section 2.2.3). In some flows the temperature and velocity gradients can be similar (the 
Reynolds analogy, Section 3.4.6), but this is not always the case as discussed further for liquid 
metals in Volume 5 (Section 2.1.2). 


In summary, the convective heat transfer at walls is strongly influenced by the detailed near-wall 
flow field, and can therefore be affected by many factors. 


Natural, Forced and Mixed Convection 


These terms describe different mechanisms of convective heat transfer that might occur depending 
on the significance of local buoyancy effects. Heat transfer by natural, forced and mixed convection 
is discussed in Volume 3 (Section 2.2.1). 


Heat Transfer Coefficients (HTCs) 


Heat transfer coefficients, or ‘film coefficients’ are extremely convenient and widely used across 
industry. Among their many advantages, HTCs: 


¢ Enable engineers to compare and discuss work on convective heat transfer, and are often 
easy to sense-check against values in textbooks or wider literature. 


« Are invaluable as an input to system performance models (Section 3.2), which may be based 
on detailed Computational Fluid Dynamics (CFD) analysis (Section 3.4), and transferring 
information about convective heat transfer between different analyses (such as from a CFD 
assessment to a structural assessment, Section 2.4). 


However, HTCs have some significant limitations, which include: 


¢ Even in relatively simple situations that involve fully developed flow over a surface, there may 
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still be significant sources of uncertainty in the HTC due to surface roughness (Section 2.2.1) 
or non-idealised boundary conditions (i.e. not isothermal or constant heat flux). 


* In practice, many real flows are complex and three-dimensional with non-idealised boundary 
conditions. This is a particular consideration and source of uncertainty/error where HTC 
correlations (Rohsenow et a/., 1998) are used in a system code (Section 3.2). 


* Engineering flows often do not have a clearly identifiable local bulk temperature (Ts,.0; Ts); 
and so skilled engineering judgement is often needed to identify, and most importantly record, 
the most appropriate approach to calculate values of HTC. For example, applying a global 
reference temperature to get a local HTC in a system where the bulk fluid temperature 
changes significantly will lead to unrealistic values of HTC. This is particularly significant 
where CFD derived flow fields are being post-processed to provide inputs to other spatial 
models, and is considered further in Section 3.4.6. 


In summary, HTCs are a useful tool. However, it is important to recognise that they are coefficients 
that attempt to represent the complex reality of heat transfer between fluids and solids. This heat 
transfer is driven by temperature gradients at walls, which are dependent on a wide range of factors, 
and are often a function of the dynamic conditions of the system and changing physical properties. 


Local and Mean Coefficients: As noted above, the HTC (h) may vary continuously over a sur- 
face in response to changes in the flow field. For most engineering analysis, the heat transfer over 
a given area of surface is of interest, and an area-averaged HTC (h) is therefore used. In situations 
where the flow changes significantly over the surface, different area-averaged HTCs may be used 
for different regions of the surface. 


Correlations: A wide range of correlations have been developed to enable HTCs to be predicted 
in various situations. Many correlations are empirical, but some are developed analytically (partic- 
ularly for laminar flow). These correlations are valid for the particular geometries, flow conditions 
and heat transfer conditions for which they were developed, and are often based on a correlation 
for the Nusselt number (Nu): 


AL ; ; ; ; 
Nu = ke a ratio of convective to conductive heat transfer in a fluid at a boundary 
f 
Where h is the HTC, L is a characteristic length scale and kr is the thermal conductivity of the 
fluid. Like h, Nu can vary continuously across a surface. Correlations may therefore be for the local 
value of Nu, or for an area-averaged value across a surface (Nu). 


Correlations for Nu or Nu are often parameterised by several key non-dimensional groups: 


1. Reynolds number (Re), a ratio of momentum forces to viscous forces; 
2. Grashof number (Gr), a ratio of buoyancy forces to viscous forces; 
3. Prandtl number (Pr), a ratio of momentum diffusivity to thermal diffusivity. 


Different correlations are appropriate for natural, forced or mixed convection, and for internal 
and external flow, since the temperature gradients in these flow fields vary significantly (Sec- 
tion 2.1.2.1), which is discussed further in Volume 3 (Section 2.2.1). In addition, correlations are 
often only valid for particular ranges of these dimensionless groups, so using correlations outside 
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of these acceptable ranges may result in errors. Checks are often included in system codes to gen- 
erate warnings where correlations are being used outside of acceptable ranges, and these should 
be considered carefully (Section 3.2). 


In laminar convection, the Nusselt number is also sensitive to the thermal boundary conditions 
and the constants that appear in the Nusselt number correlations have different values for uniform 
wall temperature (isothermal) and uniform wall heat flux thermal boundary conditions. However in 
turbulent heat convection, thermal wall boundary conditions have less influence on Nusselt number 
and in many cases the same Nusselt number correlation can be used for both types of thermal wall 
boundary conditions. 


Many industrial applications feature complex geometries and three-dimensional or transient flow 
fields with forced, natural or mixed convection, making it difficult to identify appropriate correlations 
to use (CSNI, 2015a). These factors can motivate the use of CFD, which models the temperature 
gradients and heat fluxes directly, and does not therefore use HTCs (although HTCs are commonly 
used for external boundary conditions in CFD models, see Section 3.4). 


A wide range of correlations are available in textbooks (such as Incropera et a/., 2011) for inter- 
nal/external flow, laminar/turbulent flow, and natural, forced and mixed convection with more detail 
provided in Kakag et al. (1987) and Rohsenow et al. (1998) for specific geometries. Further in- 
formation on low Pr correlations suitable for liquid metals is provided in Volume 5 (Liquid Metal 
Thermal Hydraulics). 


Correlation Inputs: Once a suitable correlation has been chosen, the variables on which it de- 
pends must be assessed. Sources of correlations often provide guidance on how to calculate 
values used in correlations, and this should be considered if provided. However, some challenges 
posed by this may include: 


Characteristic Length Scales: In many situations determining these is straightforward (e.g. where 
the internal diameter of a pipe or stream-wise length of a plate is needed), but in some cir- 
cumstances it may not be (e.g. where geometry deviates from the simple geometry on which 
the correlation is based). While guidance may be provided with the correlation, engineering 
judgement may be required to consider what values might be appropriate in the light of the 
expected flow field. For example, hydraulic diameter is often used for non-circular ducts. 


Characteristic Flow Velocity: For external flow, engineering judgement may be required to iden- 
tify a suitable free-stream velocity that can be used. For internal flow, it may be necessary 
to calculate either a ‘superficial’ velocity based on the mass flow rate or a mass-averaged 
velocity based on the flow in the region of interest. 


Fluid Material Properties: These typically vary with pressure and temperature (Volume 3, Sec- 
tion 2.3). In many situations, the local pressure variations within the flow are small, so ma- 
terial properties can be assessed at a nominal pressure for the region of interest. However, 
this is not likely to be true for temperature where there is significant heat transfer, which is 
driven by temperature gradients: 

¢ For external flow, where temperature differences cause modest differences in material 
properties a ‘film temperature’ is commonly used: Tym = (Tw + Ts,co)/2. Iterations may 
therefore be needed to update the film temperature and associated properties based 
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on the evolving wall temperature (7,,). 


+ For internal flow, the bulk fluid temperature (T, introduced above) is often used to eval- 
uate properties for fully developed flow. In some cases a film temperature based on the 
mean of the wall temperature and bulk temperature may be used: Tym = (Tw + Ts)/2. 


« Where large differences in material properties occur, a property ratio method may be 
more appropriate (Gebhart et a/., 1988). Due to the large viscosity ratio and high Pr 
in molten salts, this is covered in more detail in Volume 6 (Molten Salt Thermal Hy- 
draulics). 

The appropriate temperature to use is part of the definition of a correlation; the material 
properties are used to evaluate the non-dimensional numbers and k is used to convert Nu 
to h. Significant inaccuracies can therefore result from using the wrong definition. 


Overall, when using a Nusselt number correlation, the length scales and velocity scales used to 
calculate the values of the relevant dimensionless numbers must be the same as defined in the 
correlation (i.e. consistency with the source is important and variation in definitions does exist). 


Thermal Resistance 


The concept of thermal resistance is introduced in Section 2.1.1.1. The thermal resistance for 
convective heat transfer at a surface (R:,) follows easily from Newton’s law of cooling, using an 


average heat transfer coefficient (h) over a given surface area (A): 


- 1 1 
Q = hA(Ty — Ts.c0) = —-(Tw —Tsco) Where Riy, = —— 
(Tw — Taco) = pe (Tw — Tao) waa, 
This is for external flow; for internal flow, Ts... is replaced by T,. This thermal resistance can 
be combined with other resistances to form an overall thermal resistance (as discussed in Sec- 
tion 2.1.1). For the simple wall shown in Figure 2.1, the overall thermal resistance would be the 
sum of the thermal resistances below (neglecting thermal radiation, Section 2.1.3): 


A convective thermal resistance on the warm surface. 

A conductive thermal resistance through the solid (Section 2.1.1). 

A contact thermal resistance between the solid and the insulation, if appropriate (Section 2.1.1). 
A conductive thermal resistance through the insulation (Section 2.1.1). 


ao PpoON = 


A convective thermal resistance on the cool surface. 


In some situations resistances are used directly, such as for fouling resistances (‘fouling factors’ or 
‘dirt factors’) on the tubes of shell and tube type heat exchangers (¢, generally provided with units 
K m? W~!). Considering a cylindrical tube with fouling on the inner (¢;) and outer (¢,) surfaces, a 
thermal resistance for radial conduction through a tube of In(r./r;)/2mLk (where r, and r; are the 
outer and inner radii respectively and L is the length of the tube) and using the outer surface area 
of the tube (A, = 27r,L) as the reference area: 


1 1 G;  In(ro/ti) Co 1 
—UA,AT wh ie ge a wat 
@ we Ga ee ae) A One Ay A, 
1 1 ro fo toln(ro/F;) 1 
that — = Riz o,Ao = = i t+ Go + = 
so that [7 th, an +¢ + k ¢ fi 
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As noted in Section 2.1.1, the area (A) may be omitted from the thermal resistance (R:,) in some 
contexts, so care is needed in calculating the overall thermal transmittance (U) and the reference 
area used should be stated clearly. In multilayer systems, if there are orders of magnitude differ- 
ences between the thermal resistance of the different layers (e.g. steel and insulation) then it can 
sometimes be appropriate to approximate it as a single layer system, controlled by the layer with 
the highest thermal resistance. 


Thermal Radiation 


All surfaces continuously emit, reflect and absorb electromagnetic radiation; the intensity and spec- 
trum of wavelengths of the emitted photons depend on the temperature of the surface®. In general, 
thermal radiation exchange at surfaces depends on: 


* The angle at which the radiation arrives or leaves. 

* The wavelength of the radiation. 

¢ The material of the surface and its condition. 

¢ Fluid transparency and properties between surfaces. 

* How directional or diffuse the incoming or reflected radiation is. 


Any medium between surfaces, except a vacuum, can also ‘participate’ in thermal radiation ex- 
change (i.e. absorb, emit and scatter radiation passing through it), although in many cases this 
is minimal (as in air). Dealing with thermal radiation in full detail is complex, and so simplifying 
assumptions are often used: 


* Surfaces considered ‘diffuse’ (i.e. emission or absorption are independent of direction), ‘gray’ 
(i.e. properties are independent of radiation wavelength) and opaque (i.e. no transmission of 
thermal radiation occurs through the surface). 


* Surfaces (or parts of surfaces) considered to be isothermal. 


* The medium between surfaces considered to be non-participating. 


These simplifications allow the heat transfer between surfaces to be calculated, based on the 
Stefan-Boltzmann law for the emitted heat flux due to thermal radiation from a surface (E): 


E=ceoT* 


Where T is the absolute temperature, o is the Stefan-Boltzmann constant and « is the surface 
emissivity (Section 2.1.3.2). The absorptivity of a surface (a) defines how much incident radiation 
is absorbed (with 1 — a reflected). Kirchoff’s law states that a = ¢ at each wavelength. This means 
that highly reflective surfaces also do not emit significantly. For gray surfaces, this is simplified to 
be the case across the entire spectrum, which is a sufficiently accurate approximation for most 
applications, notably excepting semi-transparent media and where insolation (e.g. heating caused 
by solar radiation) is important. For insolation, surfaces may well not be gray due to the high tem- 
perature of the sun (hence different radiation wavelengths), so that different values of absorptivity 
and emissivity are required. 


For the temperatures in most nuclear applications, emission is almost entirely in the Infrared (IR) range of the electromag- 
netic spectrum, with few very high temperature surfaces emitting in the visible range (or ‘glowing’. 
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Aspects of participating media and non-gray interactions are discussed for molten salt in Volume 
6 (Section 2.5). A more complete description of the concepts and governing equations for ther- 
mal radiation and typical engineering approximations can be found in Incropera et al. (2011), and 
substantial details can be found in more specialist texts by Modest (2013) and Howell et a/. (2016). 


The non-linearity of thermal radiation with temperature usually requires numerical or iterative solu- 
tion. Possibly for this reason, thermal radiation is often neglected from heat transfer calculations as 
a simplification. However, this may well introduce substantial inaccuracy (even at moderate tem- 
peratures) so omitting thermal radiation from heat transfer analysis should be considered carefully. 
This is particularly true in passive cooling situations where thermal radiation can suddenly become 
significant to a point where it can provide a fundamental part of a cooling duty. 


Enclosures 


For an enclosure of NV gray surfaces, analysis is complicated by the fact that each surface absorbs 
part of the incident radiation and reflects part, and incident radiation contains components from 
multiple previous surface emissions and reflections. This leads to a set of N simultaneous equa- 
tions for the net rate that energy leaves each surface (/) of the enclosure, that are framed in terms 
of the radiosity (J;) which is the thermal radiation leaving a surface (the sum of the thermal radiation 
emitted and the reflected part of the incident irradiation). By solving for the set of radiosities, the 
heat fluxes or temperatures of surfaces can be determined (more detail is available in Incropera 
et al., 2011). 

_y And 
(Ai oy hy) 

Where F;; is the view factor between the i, j pair of surfaces, which is the fraction of the thermal 
radiation leaving surface / that is intercepted by surface /. For an N-surface enclosure, there is an 
N x N view factor matrix. Values may well be non-zero for = j, so that a surface intercepts some 
of its own thermal radiation (e.g. all concave surfaces will to some extent). View factors have the 


properties of reciprocity and summation (which are useful in determining their values) so that: 


N 
AiFij = Ay Fi and S_ Fij = 
j=1 


It is common to use the special case of the two surface enclosure, which simplifies the radiosity 


equation to: 

o(T? —T) 
Lei, . ol 4 1—€ 
e1A, — AiFiz | €2Ao 


There are a number of useful simplified cases where this can be applied, such as: 


Qi = Q1 = -Q2 = 


* A small convex object in a large enclosure like a big room or the environment. In this case 
Ai/A2 © 0 and Fi2 = 1, so Q = oe, A1(T? —T>) and the emissivity of the large surroundings 
does not need to be known. As noted previously, for incident so/ar radiation the absorptivity 
of surfaces (a) may differ significantly from the emissivity (¢). 


« Facing surfaces like large parallel plates or long concentric cylinders. Figure 2.4 shows the 
heat transfer across a gap between two solids separating hot and cold regions, such as may 
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occur in a double-walled vessel, layers of multi-foil insulation, or a contact resistance (Sec- 
tion 2.1.1.1). In parallel with thermal radiation across the gap, there will either be conduction 
through a stationary fluid, or enhancement of heat transfer by buoyancy driven convection 
(the separation distance, orientation and temperature difference between the surfaces will 
determine which). In this case Ay = Az = Aand Fy2 = 1, and the emissivity of both surfaces 


is needed: 
oA(T? — TH) 
Q= TI 1. 
—+—-1 
E1 €2 
Solid Gap _ Solid 
—-» Conduction 
( ) Conduction or 
— — Cavity Convection 
aie —~» Thermal radiation 
Figure 2.4: Thermal radiation in the gap between two solids. 
Emissivities 


The surface emissivity (€) is the ratio of radiated energy compared to that which would be emitted 
from an ideal ‘black body’ at the same temperature (E = eE,). The emissivity of a surface is highly 
dependent on the material, its surface finish (roughness) or coating and the temperature. Highly 
polished metals can have values of «€ < 0.3, or even e < 0.1, whereas some of these materi- 
als when oxidised can have values of € > 0.8. Painted surfaces, ceramics and natural materials 
like wood or paper often have e > 0.9. Collections of data for a range of materials and surface 
conditions are provided in references like Modest (2013) and Howell et al. (2016). 


The emissivity of surfaces on an as-installed component may be difficult to assess because the 
surface condition may vary through life on a plant, causing a significant source of uncertainty. This 
may need to be investigated using validation data, sensitivity studies and uncertainty quantification 
(Volume 4, Confidence and Uncertainty). Validation data also has uncertainty, however, because 
instruments like pyrometers, pyranometers and IR cameras require calibration. Surfaces can be 
painted with known emissivity paints or have small swatches of material with a known emissivity 
fixed to them, to help with taking measurements. 


View Factors 


Calculating view factors between two finite surfaces with an arbitrary orientation is a mathematically 
involved process. However, attempting this from first principles is not usually necessary. Since view 
factors are geometric parameters, either a dedicated tool with a numerical algorithm to evaluate 
them could be used, a CFD code with a Surface-to-Surface (S2S) model could be employed, or 
catalogues of pre-calculated configurations could be considered (such as Modest, 2013 or Howell 
et al., 2016). 
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Thermal Resistance 


The thermal resistance due to radiative heat transfer at a surface can be derived using the de- 
composition (7; — T}) = (Ti — T2)(T1 + T2)(T? + T?) and the enclosure equations introduced in 
Section 2.1.3.1, for example (parallel plates): 


1 
Q = —(T,—T2) where Ra=( 


i 
_ oA(Th + T2)(T? + Te) 
Rtn 


1 1 
2+i-1 


This illustrates that for thermal radiation R,, depends strongly on the surface temperatures (that 
are usually being calculated) and the complexities of multi-surface enclosures in a way that was 
not the case for conduction or convection. This makes thermal resistances less convenient to use 
when thermal radiation is involved. 


Surface Modifications 


The characteristics of the solid surface in contact with a fluid can have a significant impact on local 
heat transfer. It is possible to significantly increase or decrease the heat transfer by convection 
and/or thermal radiation by making deliberate modifications to the surface finish. It is also true, 
across industry, that surfaces rarely remain in their ‘as manufactured’ state once a component is in 
service. In order to predict the heat transfer between the solid and the fluid it is therefore important 
to consider the details and condition of the solid surface. This section discusses common surface 
characteristics and the implications for thermal hydraulic analysis. 


Surface Roughness and Heat Transfer Enhancement 


Surface roughness has the potential to significantly affect both the pressure loss and heat transfer 
through a component. A high surface roughness disrupts the viscous and thermal sublayers within 
the boundary layer. This promotes transition to turbulence and subsequent turbulence generation 
thereby increasing mixing and heat transfer. Surface roughness may be as an incidental conse- 
quence of manufacturing (e.g. machining, casting or welding) or may be deliberately controlled to 
achieve a specific heat transfer and hydraulic performance. It should also be noted that the rough- 
ness of a surface may change over time once a component is in service due to erosion or corrosion 
damage or the deposition of substances on the surface. 


The extent to which roughness affects the flow depends on the height of the roughness compared 
to the thickness of the viscous sublayer in turbulent flow (Figure 2.5). For example, the thickness of 
the viscous sublayer decreases with increasing Reynolds number and therefore a smaller surface 
roughness will have an effect. For fully developed turbulent flow in circular ducts this is illustrated 
by the well-known ‘Moody Diagram’, which shows the friction factor as a function of roughness and 
Reynolds number (Section 3.4.2.5). 


For internal flow in simple geometries, such as pipes, ducts or between plates, correlations for 
the Nusselt number under a range of conditions can be found in text books such as Rohsenow 
et al. (1998)®. Provision for random ‘sand grain’ rough surfaces is also made in most mainstream 


8 For fluids where Pr < 1, general correlations are often not applicable as discussed in Volume 5. 
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€ = Roughness height 


6, = Viscous sublayer height 


(a) Smooth Wall (b) Rough Wall (c) Complete Turbulence 


Figure 2.5: Illustration of impact of roughness on boundary layer. 


CFD codes, but careful reference should be made to the code documentation, especially when 
predicting surface to fluid heat transfer in fluids where Pr is significantly different from 1. 


The surface heat transfer from a rough surface can be several times that of a hydraulically smooth 
surface. For this reason designers sometimes deliberately roughen surfaces to increase heat trans- 
fer. Surfaces that are artificially roughened for the purpose of enhancing heat transfer often employ 
a regular pattern of protuberances or grooves (such as the ribs that may be used on fuel cladding). 
The heat transfer (and hydrodynamic) performance of these surfaces can be more complicated 
to predict. This is because the heat transfer enhancement achieved is often partially due to an 
increase in the surface area (the fin effect) as well as increased turbulent mixing in the near-wall 
flow. The performance of a specific surface design is often experimentally measured and can then 
be used as an input to predictive analysis. 


Other design modifications that may be used to enhance single-phase heat transfer include ex- 
tended surfaces (such as fins) and devices placed either on or close to the surface that change 
the flow to promote heat transfer. These are effective in both laminar and turbulent flow. Heat 
exchanger designers often make extensive use of a complex (and sometimes proprietary) combi- 
nation of these modifications to achieve a high heat transfer in a compact design. Further details 
on a greater variety of surface heat transfer enhancement techniques can be found in Rohsenow 
et al. (1998, Chapter 11) or Saha et al. (2016). 


Extended surfaces work by increasing the surface area available for heat transfer. Conjugate heat 
transfer is of vital importance in the prediction of heat transfer for a finned surface, as the temper- 
ature of the fin can vary significantly between the root and the tip. Helical fins act to both increase 
the surface area and modify the flow by introducing swirl. This increases the effective flow path and 
promotes mixing via secondary flows. 


The spacer grids used to maintain the geometry of fuel assemblies are an example of flow mod- 
ification to enhance heat transfer using a device remote from the surface of interest. These grids 
may be designed with vanes to increase flow mixing, in addition to their primary, structural support 
function. The corresponding increase in heat transfer between the fuel and the primary circuit fluid 
can be observed in the region downstream of each spacer grid. 
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Surface Deposition 


The deposition of solids on surfaces within a nuclear reactor is almost always driven by a chemical 
reaction. The chemical composition of the exposed solid surface and the fluid determine the com- 
position of the deposit and this therefore varies significantly between reactor technologies. Some 
surface deposits are beneficial, for example, stainless steel is designed to maintain a stable oxide 
layer giving it corrosion resistance. The rate of accumulation of the deposit and whether it adheres 
to the surface or breaks away are often highly dependent on local thermal hydraulic considerations, 
such as temperature and flow rate. 


All reactor technologies have the potential to form deposits in the primary circuit. That most familiar 
to UK nuclear engineers is the deposition of carbon within the primary circuit of Advanced Gas- 
cooled Reactors (AGRs) (the chemistry of which is briefly described in ONR, 2019). The carbon 
deposition occurs throughout the circuit, but especially in the higher temperature regions, including 
on the fuel cladding (Figure 2.6). The deposit reduces the effectiveness of the heat transfer from 
the fuel to the coolant in two ways. Firstly, by reducing the effectiveness of the ribs on the surface of 
fuel (by accumulating between them) and, secondly, by introducing an additional thermal resistance 
as heat from the fuel must conduct through the deposit. 


_ 


Figure 2.6: Picture of high levels of AGR fuel pin deposition (Mowforth et a/., 2015). 
This image is copyright of EDF Energy, is used with permission, and is not covered by 
the creative commons license defined in the legal statement for the present document. 


It is important for CHT modelling to consider if the presence of a deposit must be taken into account. 
This is often quite challenging as the thickness, location and thermal properties of the deposit may 
be subject to considerable uncertainty. Some deposits have complex, porous structures and some 
vary in thickness over the life of the reactor (e.g. after 40 years), building up and then breaking off 
into the flow or reducing in thickness due to a chemistry or temperature change. 


In addition to the surface-to-fluid heat transfer, changes to the surface emissivity may also need 
to be considered. Predicting the occurrence of a deposit in a quantitative manner is even more 
complex and requires a multiphysics model coupling chemistry with thermal hydraulics. Even with 
modern modelling techniques, much still needs to be determined by inspection and experiment. 


A complete study of deposits found in nuclear reactors is beyond the scope of this document. 
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However, considering the technologies most likely to comprise the next generation of NPP, the 
following sections consider relevant examples with references to enable further reading. 


CRUD 


The most important and highly studied deposit in a Pressurised Water Reactor (PWR) reactor is 
often referred to as Chalk River Unidentified Deposits (CRUD). It is used to describe the chemically 
and structurally complex scale that accumulates on the hottest regions of the reactor, including the 
fuel cladding, and is thought to be caused by sub-cooled boiling that precipitates solutes in the pri- 
mary circuit fluids on the cladding surface. At typical PWR temperatures this process is aggravated 
by the inverse relationship between the solubility of key compounds and the water temperature. The 
CRUD accumulates as a porous layer, with boiling (and therefore further precipitation) occurring 
within its structure. The resulting CRUD is chemically complex, comprising a variety of compounds 
including both primary circuit constituents and corrosion products. A similar, although chemically 
slightly different deposit is also observed in Boiling Water Reactors (BWRs). 


CRUD is of concern for several reasons including its effect on the core neutronics. Compounds 
of boron (intentionally present in the primary circuit water as boric acid used to control reactivity) 
can accumulate in the CRUD, locally suppressing the neutron flux. This can cause a phenomenon 
called CRUD Induced Power Shift (CIPS), also called the Axial Off-set Anomaly (AOA), where a 
downward axial shift in power distribution is observed. Control rods enter the core of a PWR from 
the top, so this could cause a delay in reactivity control. By the time the fuel is removed from the 
reactor, the CRUD itself can contain activated isotopes leading to increased radioactivity of the 
fuel. From a thermal hydraulic perspective, CRUD may impair the heat transfer from the fuel to the 
primary circuit fluid. The locally increased fuel cladding temperatures can lead to an increased rate 
of corrosion damage to the cladding. The presence of CRUD can also increases the fuel rod friction 
coefficient, increasing the pressure drop in affected areas of the core and potentially altering the 
core flow distribution. 


Due to its importance to current generation reactor performance, CRUD has been extensively stud- 
ied and significant guidance is captured in documents such as EPRI (2004), EPRI (2014) and EPRI 
(2020). Thermal hydraulics is of significance both to the prediction of the formation of CRUD and 
its impact on heat transfer. The prediction of CRUD formation requires the coupling of chemistry 
and thermal hydraulics (as described in Henshaw et a/., 2006). Multiphysics nuclear codes have 
been developed for this purpose, for example, the MAMBA code as described in Kendrick et al. 
(2013) was incorporated with a whole core model to predict CIPS using VERA, as part of the 
CASL programme (Collins et al., 2018). 


Prediction of heat transfer through a CRUD layer is complicated by the boiling that can occur within 
the porous structure and this is further complicated by the effect of the concentration of solutes 
on the saturation temperature of water. EPRI (2011) describes experimental work to determine 
the effective thermal conductivity of CRUD carried out using the Westinghouse Advanced Loop 
Tester (WALT) facility. The effective thermal conductivity is found to be highly dependent on the 
heat flux in addition to the thickness of the deposit. With increasing heat flux the fluid within the 
CRUD undergoes single-phase heat transfer, nucleate and ‘wick’ boiling, and dry-out. The study 
also found that under specific conditions of a small (around 15 um) CRUD thickness and moderate 
heat flux, the CRUD acted to enhance heat transfer from the surface when compared to a clean rod, 


18 of 84 


Volume 2 


2.2.2.2 


2.2.2.3 


Technical Context 


by promoting the initiation of sub-cooled nucleate boiling. However, at higher CRUD thicknesses or 
higher (or lower) heat fluxes the CRUD is generally found to impair the heat transfer. 


Oxide Layers in Heavy Metal Reactors 


A protective layer on the surface of steel reactor components is required in a heavy liquid metal 
environment. Without it, the steel is attacked by dissolution of its components into the lead or 
Lead-Bismuth Eutectic (LBE) primary coolant. A popular solution to this problem is to introduce 
oxygen into the coolant to promote and maintain the formation of thin oxide layers on the steel 
surfaces. The level of oxygen needs to be carefully controlled: too little and the oxide layer will not 
be maintained; too high and the oxide layer may grow to a size where it spalls off the surface or 
oxides may form in the coolant itself (e.g. lead oxide). Excess oxide may then be transported by 
the coolant and deposit in the cooler regions of the reactor, potentially disrupting flow paths. 


The nature of the oxide layer depends on the constituents of the steel. Oxides of iron and chromium 
are often present and a complex two or three layer structure can develop (further described in NSC, 
2015). Due to the importance of corrosion resistance in heavy liquid metal reactors, an extensive 
range of relevant experiments have been published using a wide range of steels (for example 
Deloffre et a/., 2002). However, it is noted that work is still ongoing with no single best way forward 
identified to date. 


The effect of the oxide layer on the conjugate heat transfer depends on its thickness and conduc- 
tivity. The thermal conductivity of the oxide layer will depend on its composition and structure (e.g. 
the proportion of voids within the oxide). It is extremely likely to be lower than that of either the 
steel or the coolant and will therefore reduce the surface-to-fluid heat transfer, however, the effect 
may be negligible if the layer is thin. Despite the amount of work that has been undertaken, oxide 
layer studies focus on the formation, composition and structural characteristics of the layer (i.e. 
the factors that matter most to corrosion resistance) rather than the thermal properties. There is, 
therefore, likely to be considerable uncertainty in the thermal properties of the layer. 


Metallic Deposits in Sodium and Molten Salt Cooled Reactors 


The surface deposition within Sodium-cooled Fast Reactors (SFRs) or Molten Salt Reactors (MSRs) 
differs from other technologies. Of particular significance is that the oxide layer normally present 
on the surfaces of structural steels is not stable in these environments. The bare metal is therefore 
exposed to the coolant, potentially making it more vulnerable to corrosive attack. 


Under these circumstances, alloying agents within the steel, such as chromium or nickel, are dis- 
solved into the coolant in hotter regions of the reactor and these can then be deposited (or ‘plate 
out’) in cooler regions, such as heat exchangers. Any contaminants present in the coolant (in- 
evitable in low concentrations) can also lead to the formation of oxide, carbide or nitride deposits. 
For MSR designs where the fuel and therefore the fission products are dissolved within the primary 
circuit fluid, deposition of some fission products has been found to occur on reactor surfaces. 


Surface deposition is not an area of high concern for these reactor technologies. However, over 
time, the level of deposition may result in flow restrictions in narrow passages or heat transfer 
impairment (Pioro, 2016). 
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Graphite Dust in High Temperature Gas Reactors 


Where graphite is used as a moderator, there is the potential for the generation, dispersion and 
deposition of this dust on surfaces within the reactor. This is of particular concern for pebble bed re- 
actors where the movement of the pebbles generates dust by abrasion. This dust is then deposited 
on the fuel itself or may be transported by the primary circuit fluid and deposited on other reactor 
components (for example, the steam generator, as described in Peng et al., 2016). The dust is 
radioactive so, as well as potentially impeding heat transfer, it complicates maintenance activities 
and could contribute to a radiological release in the event of a breach of the primary circuit. 


Surface Corrosion and Erosion 


The International Organization for Standardization (ISO) definition of corrosion is ‘a physiochemical 
interaction leading to a significant deterioration of the functional properties of a material, or the 
environment with which it has interacted, or both of these’ (Féron, 2012). Corrosion damage can 
be found in a variety of engineering environments and applications, with one of the most widely 
encountered being the corrosion of metals in the presence of water. 


Erosion is the degradation of a surface due to mechanical action. In fluid mechanics this can be by 
the impingement of a fluid on a surface or abrasion by solid particles, gas bubbles or liquid droplets 
suspended in the fluid. Corrosion and erosion can change the roughness of a surface and/or the 
local material properties of the surface material. Therefore, the possibility of corrosion and erosion 
damage to a surface needs to be taken into account when considering heat transfer. 


Thermal hydraulics is also important in the occurrence of corrosion and erosion. For erosion, the 
momentum of the flow (and any particles suspended in it) is directly responsible for the damage. 
Important factors in the prediction of the rate of erosion include the flow velocity and the angle 
of impact with the surface. For corrosion, the fluid flow may transport the reactants involved to 
the surface and the rate of the chemical reactions associated with corrosion is often tempera- 
ture dependent. Erosion and corrosion can also act together to cause more damage that either 
would alone. For example, in Flow-Accelerated Corrosion (FAC), fast flowing fluid is responsible for 
stripping away (or dissolving) corrosion products from the surface leaving it open to further attack 
(Figure 2.7). 


Fluid o——— 


Metal 


(a) Diffusion (b) Erosion (c) Dissolution 
Figure 2.7: Illustration of processes occurring during Flow-Accelerated Corrosion. 


In NPPs, the control of damage to reactor components is of key importance to both nuclear safety 
and the economic performance of the plant. Difficulties in the replacement of components in a 
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nuclear environment mean that the selection of materials with a long (>40 years) lifespan in the 
relevant (and industrially unique) environments are required. Additionally, irradiation can act to 
change the properties and structure of materials over their lifetime making them more vulnerable 
to damage. 


Corrosion and erosion damage to NPP components is a large topic involving material science and 
chemistry as well as thermal hydraulics. Féron (2012) provides a useful overview of the subject of 
corrosion, including corrosion mechanisms commonly found in NPPs and in-service inspection an 
monitoring. The remainder of this section considers only aspects that potentially impact on CHT 
modelling rather than the wider topic. 


Surface Geometry 


An important consideration when attempting to model a corroded surface is that some corrosion 
products have a different physical structure and volume to the original metal material, as well as 
different material properties. The changes to the geometry of the corroded component may have 
a significant affect on the fluid to surface heat transfer if, for example, a geometrical heat transfer 
enhancement feature (see Section 2.2.1) has been destroyed or covered by oxide. 


Aqueous corrosion is a commonly found example where the corrosion products (e.g. aqueous iron 
oxide) can have a significantly different structure and volume to the original steel. The majority of 
NPPs in current operation are Light Water Reactors (LWRs) that typically use water in both the 
primary and secondary circuits. These are a mature technology, so have had the advantage of 
learning from the reactors of the past with regard to material selection, but potentially higher risk 
areas are still present. The biggest corrosion issues are often found on the secondary (steam side) 
components (Féron, 2012), including the primary heat exchanger, and are often associated with 
FAC. This is likely to be a concern for all reactor technologies with a steam secondary system. 
Additionally, peripheral reactor systems often include pipes and tanks that are largely unprotected 
from the weather and so are vulnerable to damage by corrosion. 


Surface Roughness 


One of the easiest surface changes to envisage due to damage by erosion and corrosion is the 
surface roughness. For a smooth surface, damage can roughen it, enhancing heat transfer and 
increasing the friction factor (potentially resulting in higher pressure losses) as described in Sec- 
tion 2.2.1. It is also noted that surface roughness can have a very significant effect on initiation and 
progression of corrosion, with rougher surfaces often being more vulnerable that smoother ones. 
For example, Platt et a/. (2015) illustrates the effect of surface roughness on oxidation of the fuel 
cladding material Zircaloy-4. 


Corrosion and erosion can also act to reduce the surface roughness where ribs or other heat 
transfer enhancing features are present. This has significant overlap with the consequences of 
deposition as discussed in Section 2.2.2. 
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Surface Emissivity 


Changes to a surface finish have the potential to significantly affect the radiative heat transfer by 
changing the surface emissivity. For example, polished stainless steel has an emissivity of around 
0.3 whereas oxidised stainless steel has an emissivity of closer to 0.8 (Incropera et al., 2011). As 
described in Section 2.1.3, the transfer of energy by thermal radiation is proportional to T*, so 
changes to emissivity can significantly impact radiative heat transfer in nuclear reactors where the 
temperatures are high. This is important for passive cooling systems that rely on thermal radiation, 
as surface finish degradation could reduce the heat transfer effectiveness. 


For commonly used materials, such as steel, values of emissivity are often readily available for a 
variety of surface conditions. However, a full range of information may be harder to find where spe- 
cialist materials are used or where erosion or corrosion has degraded a surface coating, exposing 
the material underneath. 


Surface Material Properties 


Corrosion in particular can significantly alter the thermal and structural properties of the surface 
material. Generally, corrosion products (including oxide layers) are not considered to contribute to 
the structural properties of the component, so corrosion of a pipe (say) is effectively a thinning of 
its walls. Corrosion of alloys can alter the material properties of the surface without producing a 
deposit. 


Fluids such as lead, LBE and molten salts can leach alloying agents such as nickel and chromium 
from steels, especially along the steel grain boundaries. This leaves a layer of material with a 
weakened, more porous structure close to the surface. Further information on steel corrosion in 
molten salt and heavy liquid metal environments is given in Zheng and Sridharan (2018) for salts 
and Zhang (2009) for lead/LBE. 


For CHT, the specific material properties (especially conductivity) of any corrosion products present 
onthe surface may be important. In a similar manner to surface deposition, the specific composition 
and material structure of any corrosion products may be hard to determine/predict with a high level 
of certainty, and so the information needed to determine thermal conduction through them may not 
be easy to acquire. 


Surface Coatings 


In engineering applications it is not unusual to deliberately apply coatings to the surface of compo- 
nents, paint being a very simple example. These coatings are often applied to increase the resis- 
tance of the surface to damage, but can be used for other reasons (e.g. to change the appearance 
of the surface). 


A relevant recent example of coatings considered for a NPP application is in the development of 
accident tolerant fuels for LWRs as described in Terrani (2018). These are intended to reduce the 
rate of oxidation of the underlying Zr fuel cladding in the presence of high pressure steam under 
severe accident conditions. It is specifically noted in Terrani (2018) that, in this case, the coating is 
considered sufficiently thin to have no significant effect on the heat transfer. 
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An advantage of deliberately applied coatings (as opposed to surface deposits or corrosion) is 
that their effect on the thermal performance of a component can be considered (and potentially 
measured) at the outset. The information required to take account of these coatings in any CHT 
modelling is, therefore, more likely to be available. 


Material Properties 


Temperature gradients within a solid and how quickly the temperature changes over time is de- 
pendent on the thermophysical properties of the solid (density, specific heat capacity, thermal con- 
ductivity and emissivity). For transient simulations all of these properties are important in order to 
account for the thermal mass of the solid. In general, the properties of solids vary due to a range 
of factors that should be considered, including: 


Nature of Solid: The properties of solids may be altered by a range of factors (such as contam- 
ination, addition of different species and ageing due to radioactivity). This means that the 
material properties may change through the life of the plant, particularly for those compo- 
nents in high radiation areas, such as the core. The key factors to consider and their impact 
are likely to be case specific. 


Temperature: Temperature changes often cause large variations in solid properties. Therefore, 
the variations in properties with temperature are a key aspect to consider when performing 
thermal analyses with conjugate heat transfer. As for fluid properties, temperature varying 
properties can be implemented by linearly interpolating or fitting measured data (Volume 3, 
Section 2.3.2). 


Surface Finish: As discussed in Section 2.1.3.2, the emissivity of a surface is highly dependent 
on the material, its surface finish and the temperature. Therefore, the emissivity of a solid 
can change through the life of the plant as the surface becomes oxidised or covered by a thin 
layer of deposit. As a result, the emissivity is a significant source of uncertainty in material 
property data that may need to be investigated. 


Porosity: Some materials are porous (contain pores/voids that are filled with fluid), such as in- 
sulation or graphite. The impact of the fluid on the overall material properties needs to be 
considered and taken into account, if necessary. However, this is essential for porous zones 
(regions used to represent solid and fluid components, Section 3.3.4), such as perforated 
plates or tube banks, and effective material properties for the porous zone need to be care- 
fully calculated (Section 3.4.3.3). 


Anisotropy: Some materials, such as composites, reinforced concrete and porous media (Sec- 
tion 3.4.3.3), also have anisotropic properties, such as thermal conductivity, which may 
change over time. This can significantly impact the thermal gradients within the solid and 
so needs to be carefully considered. Anisotropic thermal conductivity is available in most 
CFD codes, but requires knowledge of the material specification or can be calculated using 
a detailed thermal model of the material. 


The rest of this section identifies some sources of material property data and considers how these 
might be used. However, it is also important to understand the accuracy of the data itself, and so 
whether this needs to be considered as part of a sensitivity analysis or uncertainty quantification. 
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The properties of fluids are discussed in Volume 3 (Section 2.3). 


Sources of Property Data 


Possible sources of material property data for a number of common solids are presented below, 
although general values can be found in most text books, such as Incropera et al., 2011. Whatever 
source is chosen, this is a key input to analysis work, so the accuracy of the data should be 
considered and the approach used documented. 


Steel: Most SSCs are manufactured from different types (carbon, alloy or stainless) and grades 
(classification based on alloy content and manufacturing process) of steel. The type/grade of steel 
can significantly affect the thermophysical properties of the material, and so it is important to iden- 
tify the actual steel used and determine its properties. It is worth noting that the thermal conductivity 
of stainless steel is significantly lower than carbon steel. Published sources of data for common 
steels include ASME (2019a), text books such as Davis (1998) or alternatively specification sheets 
may be obtained from the steel manufacturer, if known. In addition, IAEA (2006) and IAEA (2009) 
provide properties of some steels, as well as fuel material, cladding and moderator properties. 


Insulation: Insulation is used inside and outside containment to reduce heat loss from hot sur- 
faces. This is designed to have a low thermal conductivity and covers a wide range of different 
materials, such as Reflective Metallic Insulation (RMI), encapsulated fibreglass and calcium sili- 
cate. Published sources of data for common insulation materials include Incropera et al. (2011), 
alternatively specification sheets may be obtained from the manufacturer, if known. The uncer- 
tainty in insulation properties can be significant due to poor fitting, compression, ageing/damage 
over time, gaps between sheets or penetrations through the insulation, and so the effective thermal 
conductivity may be much higher than expected. 


Concrete: Concrete is the most commonly used civil engineering construction material in the 
nuclear industry with different mixes used depending on the situation and application (density, 
strength and durability), and may be reinforced with steel rebar or internal cables (prestressed). 
The properties depend on the concrete mix used and can change over time due to ageing and 
thermal effects. Published sources of concrete material data include ORNL (2006) and The Con- 
crete Centre’s dynamic thermal properties calculator’. 


Structural Integrity 


As discussed in Section 1.2.2, CHT analysis is often used to provide temperature fields in support 
of structural integrity assessments. As well as the direct impact of spatially varying temperature 
fields on stresses in materials, temperature fields are also critical in defining the material proper- 
ties that are used in assessments. The properties of metals (such as Young’s Modulus, yield stress 
and fracture toughness) change with temperature, and therefore a good resolution of the temper- 
ature field is key to a well-constructed and rigorous structural integrity assessment. This section 
introduces some key structural integrity concepts. 


www.concretecentre.com/Publications-Software/Design-tools-and-software/Dynamic- Thermal-Properties-Calculator- (1) 
-ASPX 


24 of 84 


2.4.1 


2.4.2 


Technical Context 


Loading Categories 


Before discussing how temperature fields are used in structural integrity assessments, it is neces- 
sary to introduce some concepts with respect to loading and consequently stress. Loading and the 
resulting stresses are categorised as either primary or secondary: 


Primary loads are stresses that are required to equilibrate applied loads. They are known as 
‘sustained’ or ‘persistent’. As a structure deforms, primary loads do not relax (diminish) or 
reduce their effect on the structure, and may therefore cause unbounded deformation leading 
to failure. Examples of primary loading include loads arising from gravity or internal pressure. 


Secondary loads are stresses that diminish as a structure deforms. An example of secondary 
loading is a load arising from constrained thermal expansion. The stresses are generated 
by the constrained expansion, so if the constraining structure were to displace, the stresses 
within the constrained structure would typically relax and not drive the structure to failure 
alone (although there may be complex interactions with other mechanisms). 


While the above example categorises thermal stresses as secondary, in practice thermal stresses 
can be either secondary or ‘long range’ (primary in how they are assessed). A through-wall tem- 
perature gradient in a pipe and the resulting bending stress would typically be categorised as 
secondary. However, a long length of expanding pipework would be categorised as long range, 
when considering its effect on a connecting nozzle, anchor or vessel. This is because the defor- 
mation of the nozzle, anchor or vessel, will likely not lead to sufficient relaxation of this thermal 
expansion to avoid driving failure in the component. 


Assessment Concepts 


The purpose of a structural integrity assessment is to limit the possibility of failure to within recog- 
nised acceptable limits. The assessment procedures are numerously codified and these codes, 
although based on the same principles, are defined by national bodies and often have a basis in 
legislation. Examples of design codes used to design of piping and pressure vessels include ASME 
(2019b), BS EN (BSI, 2014, BSI, 2017) or AFCEN (AFCEN, 2018a, AFCEN, 2018b). These codes 
relate to new design, but the codified methods are used in assessing existing components as well. 


Failure of an SSC would typically be defined as a breach of a pressure boundary or an inability 
to perform its intended function, both of which may have significant safety and economic impact. 
Protection against failure from the following broad concepts are assessed: 

1. Plastic collapse: failure through gross plastic deformation. 


2. Gross plastic strain: large stresses that result in plastic strains that are greater than the 
material ductility. 


3. Cyclic loading: alternating loads that typically result in fatigue and/or alternating plasticity 
resulting in a ‘ratchet’ failure. 


4. Buckling: instability in structures due to compressive loading effects. 


5. Creep: this includes mechanisms like creep rupture, creep fatigue or creep crack growth 
resulting from high temperature service. 
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6. Defect tolerance: only routinely deployed for highly safety significant structures. This con- 
siders the tolerance of a structure to known or postulated defects (or cracks) and includes 
assessment of fatigue crack growth. 


Assessment to these concepts is either conducted by rule and/or Finite Element Analysis (FEA) 
supported by assessment. The route is usually set by the assessment code, but where loading or 
geometry becomes complex (particularly with temperature fields) then the assessment is likely to 
be supported by FEA. An example FEA model is shown in Figure 2.8. 


Stress Field 


Temperature Field 


Figure 2.8: Example of stresses predicted by an FEA model. 


Temperature fields are used extensively in these assessment concepts. Of particular interest are 
complex spatially varying fields (such as Figure 2.8) in a structure resulting in differential thermal 
expansion, through-wall temperature profiles (typically from thermal shocks) that result in bending 
stress profiles, or whole system temperature gradients or differentials resulting in long range ther- 
mal stresses. Temperature gradients are a key concern, as they set up differential expansions of the 
structure and this generates stresses. However, isothermal states are of equal concern, as creep 
in high temperature isothermal states and transients between two isothermal states are important 
(e.g. differential expansion between hot and cold plant states). An overview of how temperature 
fields might be used in each of the six failure concepts is discussed below. 


Plastic Collapse: The load arising from temperature fields are only considered in plastic collapse 
assessments if they are judged to be long range. An example was given above of a long length 
of expanding pipework which could result in a primary load on any connected nozzle, anchor or 
vessel. The reason truly secondary loading is not assessed is that through yielding of the structure, 
the secondary load would relax and not lead to plastic collapse. 


26 of 84 


Volume 2 


Technical Context 


Gross Plastic Strain: The combined primary and secondary stresses are considered across 
any cycles that the component may see. The most severe of these tend to be between hot and 
cold plant states or thermal shock events. This will often see large reversing stresses that should 
be checked to ensure their total stress range does not exceed the ductility of the material of con- 
struction in one single cycle. The exhaustion of ductility leads to the rupture of the component. 


Cyclic Loading: This considers loading that is applied to structures that could, through repeated 
alternating application, drive a failure of the structure. Codified methods often have amalgamated 
methods that deal with both fatigue and ratchet in one assessment but they are two different con- 
cepts. Both of these concepts are impacted by both primary and secondary loading: 


Fatigue is the initiation of crack-like defects in the structure through the repeated application 
of loads and is typically split into high and low cycle fatigue (Figure 2.9). The boundary 
between low cycle and high cycle fatigue is typically 10,000 cycles. Above this is high cycle 
fatigue, where small elastic alternating stresses initiate a crack like defect. Below a certain 
alternating stress level (that is material, geometrically and environmentally dependent) no 
defect will be initiated no matter how many cycles occur, and this is the endurance limit. 
Low cycle fatigue is typically below 10,000 cycles and is for larger stress ranges that see 
alternating plastic stresses. Low cycle fatigue is again concerned with initiating a defect, but 
through far fewer cycles compared to high cycle fatigue. 


Ratcheting is incremental increases in plasticity brought about by alternating loads that could 
ultimately lead to ductility exhaustion of the material. If the loading is shown to be within 
the ratchet boundary then the loading is said to have ‘adapted’ and any repeat application 
of the loading will not result in any further net accumulation of plastic strain. If the repeat 
loading is not within the ratchet boundary, the repeat application of load results in incremental 
accumulation of plastic strain that ultimately leads to ductility exhaustion. 


Both of these concepts can be driven by thermal loads. Start up and shutdown of plants can be 
significant contributors, as can local thermal effects where rapid heating and cooling of the surface 
due to eddies in the adjacent fluid result in through-wall bending stresses. 


Stress Range 


Endurance Limit 


Low Cycle 
Fatigue 


High Cycle Fatigue 


1 ~100 ~10,000 Number of Cycles to Failure 


Figure 2.9: Cyclic loading of varying stress and failure mechanisms. 
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Buckling: The stability of vessel shells or piping networks is important, and buckling can be 
driven by constrained thermal expansion. The constrained thermal expansion can cause a com- 
pressive stress field that could result in buckling. 


Creep: When metals operate for long periods of time at temperatures above approximately one 
third of their melting temperature, irrecoverable strains accumulate over time. These temperature 
and time dependent strains are referred to as creep strains. Two general types of failure can occur 
due to the accumulation of creep strain: 


¢ Creep rupture results in the failure of a gross section of a component due to high accumula- 
tion of creep strain caused by the presence of primary loading. 


* Creep ductility exhaustion can lead to localised cracking and is generally caused by creep 
strain accumulation from the relaxation of secondary stresses cause by thermal fields. 


Where both creep strains and cyclic loading occur in a structure they can become synergistic, a 
phenomenon referred to as creep-fatigue. The presence of both simultaneously can lead to crack- 
ing of the component in unexpectedly short timescales than if fatigue and creep were assessed 
separately. Alternating thermal stresses, or creep strain from the relaxation of thermal stresses can 
have a significantly deleterious effect on creep-fatigue life. 


If a structure operating in the creep temperature range is cracked, the crack will grow in the absence 
of any fatigue loading due to the accumulation of creep damage at the crack tip. If fatigue loading 
is also a contributor to crack growth, synergy between the fatigue and creep mechanisms can lead 
to accelerated crack growth rates. Alternating thermal stresses, or creep strain from the relaxation 
of thermal stresses can be a significant contributing factor in creep crack growth rates. 


Defect Tolerance: In highly safety classified components there is likely to be an assessment 
of known or postulated defects. Any metallic structure can have a defect present at manufacture. 
Components are inspected using Non-Destructive Testing (NDT) techniques when they are built, 
but they typically have a threshold, below which defects cannot be detected. This threshold is driven 
by the physical constraints of the method deployed and the skill of the operator. In high safety duty 
components, these inspections are repeated at certain intervals through their life to identify defects 
that may have appeared through fatigue or environmental degradation, or even identify an original 
manufacturing defect that was missed. 


Defect tolerance seeks to assess the tolerance of a structure to the presence of known or postu- 
lated defects. Defect tolerance assessments are based on fracture mechanics principles and de- 
termine the proximity of any defect to failure through either fracture or plastic collapse. A detailed 
appreciation of the stress field around these defects is key to these assessments and thermal fields 
are no less important. The loads will be classified as primary or secondary, but their impact can be 
significant in either case. 


In situations with complex thermal fields and little margin, FEA can be deployed to model a crack 
sitting within a thermal field. Key crack parameters (J integrals) would then be evaluated and the 
impact of the thermal field assessed. 


The growth of any defect can also be considered through an assessment of fatigue crack growth. 
This is again driven by a detailed understanding of alternating stresses within the structure and 
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thermal fields resulting from differential thermal expansion, through-wall effects and long range 
effects would be considered. 


CHT Applications 


This section introduces some examples of applications of using CHT analysis to support structural 
integrity assessments, while the impact of buoyancy on surface temperatures is discussed in Vol- 
ume 3 (Section 2.1.3). Particular focus is given to understanding temperature fluctuations in solids 
that could lead to thermal fatigue. 


Thermal fatigue: High-cycle thermal fatigue in T-junctions was initially investigated in the context 
of Liquid Metal-cooled Fast Reactors (LMFRs) in the 1980s, but received increased interest (Kuhn 
et al., 2010, Selvam et al., 2015 and Tunstall et a/., 2016) following the incident at the Civaux-1 plant 
in France. The Committee on the Safety of Nuclear Installations (CSNI) Vattenfall T-Junction bench- 
mark (CSNI, 2011) investigated the ability of CFD to predict the amplitude and frequency of the 
thermal fluctuations, and demonstrated the improved performance of Large Eddy Simulation (LES) 
to simulate this phenomenon compared to Unsteady Reynolds-Averaged Navier-Stokes (URANS). 
CSNI (2015b) gives details on the frequencies of interest for high-frequency thermal fluctuation (or 
‘thermal striping’), where failures have been seen and general best practice guidance. 


Pressurised Thermal Shock (PTS): Thermal loads on an RPV in a PWR may lead to very large 
stresses and growth of any flaws present. Cold water injection into the downcomer of the RPV can 
occur during Emergency Core Cooling System (ECCS) operation in response to a Loss-Of-Coolant 
Accident (LOCA). Therefore, PTS is a key reactor safety issue that may impact the lifespan (and 
subsequent economics) of the plant due to radiation embrittlement of materials over time (IAEA, 
2010). This has been investigated extensively, and has been simulated in the most recent CSNI 
CFD benchmark on cold leg mixing’. The process of performing and validating CFD simulations 
relevant to PTS is also described in Study D: System Code and CFD Analysis for a Light Water 
Small Modular Reactor. 


Thermo-mechanical loads: Thermo-mechanical studies of components are undertaken to cal- 
culate solid temperatures in order to determine the mechanical loads. One such example is a 
thermo-mechanical study of the bolts holding the PWR core shielding (Rupp et a/., 2009). The 
objective was the evaluation of the mechanical integrity of 1,000 bolts used to hold the peripheral 
thermal shielding part of the nuclear core, through a structural mechanics analysis of the ther- 
mal stresses. A coupled CFD/heat conduction simulation, using two separate codes, one for CFD 
(Code_Saturne) and one for heat conduction (SYRTHES) was carried out. The simulated tem- 
perature field of the bolts holding the peripheral shielding was transferred to a structural analysis 
code. 


8 The report for which has not been released at time of writing. 
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Fuel cladding temperature: The maximum allowable temperature of fuel cladding under steady- 
state, transient and accident conditions can limit reactor performance (IAEA, 2008). This is appli- 
cable to standard and wire-wrapped fuel assemblies, and is closely coupled to the core neutronics. 
This is a key output from subchannel code simulations of the core behaviour, but has more recently 
been investigated using CFD analysis including benchmark simulations of a 19-pin wire-wrapped 
rod bundle (Merzari et a/., 2016). 


Modelling Challenges 


This section introduces the typical challenges that may be encountered in modelling heat transfer 
in solids and fluids. Addressing these challenges is considered in Section 3. 


Timescales: There is generally a large disparity between the timescales associated with detailed 
flow effects (less than a second), solid conduction (hours) and the whole system (multiple hours or 
days). Although a system or subchannel code can simulate long transient scenarios, they do not 
capture all complex flow features. By comparison, a CFD model can resolve the flow effects, but 
the large disparity between the characteristic timescales of the fluid flow (e.g. vortex shedding) and 
solid conduction (e.g. thermal mass) means that transient modelling of a CHT problem is impracti- 
cal for large components or long transients, without a carefully designed modelling approach. 


Heat Transfer Coefficients: HTCs, as discussed in Section 2.1.2.3, are widely used across 
industry as inputs for system codes, external boundary conditions for CFD models and to transfer 
information between models. However, empirical correlations are only applicable to simple flows 
with idealised boundary conditions, and it is often difficult to extract appropriate values of HTC 
from a CFD simulation. Therefore, there can be significant uncertainty in both empirical and CFD 
calculated values of HTC. 


Contact Resistance: It is often assumed for analysis that two solid surfaces are in perfect con- 
tact with each other, but as discussed in Section 2.1.1.1 this is rarely the case. The imperfections 
and roughness on solid surfaces means that the amount of contact is unknown, varies between 
components and can change with time. Therefore, this represents uncertainty in the thermal resis- 
tance associated with the conductive heat transfer that should be considered. 


Thermal Expansion: As solids heat up, the material expands based on the thermal expansion 
coefficient (1m of stainless steel expands almost 10mm between room temperature and 600 °C). 
This can potentially change the size of gaps between solids at different temperatures, e.g. the 
gap between the reactor and guard vessel for a Passive Decay Heat Removal System (PHRS). 
Since manufacturing drawings represent the cold geometry, it may be necessary but challenging 
to accurately predict the hot geometry and the impact/importance of thermal expansion needs to 
be carefully considered. 
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External Boundary Conditions: The heat lost to the environment (e.g. ambient air) can be 
sensitive to the applied external boundary conditions, as it often has a high thermal resistance. 
However, it is not possible to include an explicit representation of the external environment in 
most models, and so empirical correlations are used. This usually involves assuming natural or 
forced convection correlations based on the local component parameters. Therefore, it is often 
hard to take into account additional factors, such as varying wind conditions or buoyancy driven 
flow within an enclosure causing forced convection on a surface instead of locally calculated natural 
convection. This means that there can be significant uncertainty in empirically calculated external 
boundary conditions. 


Cluttered/Complex Geometry: A real reactor geometry is large and complex with multiple bolts, 
flanges and cluttered geometry, such as heat exchangers. In most cases, it is not possible to 
model all of the geometry in detail, and so decisions need to be made on an appropriate level 
of simplification. This will depend on the objective of the simulation, and needs to be reviewed 
within the context of a graded approach. In some cases, small features, such as fillets, can have 
a large impact on the thermal mixing at a junction. Therefore, the level of simplification needs to 
reduce around the region of interest. Models can be simplified by reducing geometrical complexity 
or representing cluttered regions as a porous zone, although the effective material properties and 
pressure drop/turbulence are often hard to calculate with confidence (Section 3.3.4). 


Heat Transfer in Turbulent Flows: Unlike in pure aerodynamics, which is dominated by the 
pressure-inertia balance, while turbulent effects are weaker, temperatures depend only on convec- 
tion and turbulent heat flux balance and so the predictions are much more sensitive first to the kine- 
matic turbulence model chosen and second to how it is extended to heat fluxes. Therefore, large 
differences can be observed in temperatures between simulations as a result of model choices. 
This means that predicting heat transfer with confidence in complex flows, such as separation and 
impingement flows, may be particularly challenging. 
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This section discusses the methods used for predicting CHT, providing a concise overview to as- 
sist engineers performing analysis. Typical coupling methods are considered, before assessments 
using system and subchannel codes, CHT approaches, CFD analysis and experimental methods 
are discussed in more detail. 


The overall approach to Nuclear Thermal Hydraulics (NTH) analysis is introduced and discussed 
in Volume 1 (Section 4.1). In general: 


* System codes are often used to predict the performance of whole systems or groups of 
systems, while CFD is generally used to understand the detailed flow field in specific com- 
ponents. 


* Application-specific codes (like subchannel or containment codes) are often used where 
some flow field modelling is needed, but a full CFD approach is too cumbersome or compu- 
tationally expensive. 


* Experiments are often used to understand overall behaviour using scaled systems, investi- 
gate detailed flow effects and provide comparative data to validate analysis methods. 


The coupling together of different methods (i.e. multiscale and multiphysics analyses) are dis- 
cussed in Section 3.1, while the different CHT approaches are considered in Section 3.3 as these 
are more relevant to CFD analysis. 


Coupling Methods 


Multiscale and multiphysics coupling is being increasingly used as part of the reactor design pro- 
cess, and offers significant advantages by combining the strengths and benefits of each analysis 
methodology. This is an active area of research across the nuclear industry to improve the integra- 
tion and coupling between codes (Section 4). 


It is important to note that Verification and Validation (V&V) should be undertaken for the coupled 
code as a whole, as well as the component models within the coupled code. 


Multiscale 


Multiscale modelling couples together multiple models at different scales and levels of resolution to 
simultaneously simulate a system. In NTH analysis, this means coupling together system (whole 
plant), subchannel (core behaviour) and CFD (detailed prediction of specific areas) codes. Zhang 
(2020) demonstrates that coupling has been achieved across most of the system/CFD codes de- 
scribed in Volume 1. 
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The most common multiscale modelling is the coupling together of system and CFD codes in order 
to simulate the behaviour of a plant/system under normal operation or fault scenarios with CFD 
models to resolve complex aspects of the flow in specific areas. This overcomes some of the limi- 
tations of system codes to model detailed mixing/stratification in pools, plena and complex geom- 
etry. In some cases, the increased computational expense associated with coupling is necessary, 
as otherwise it can lead to inaccurate system code results (Anderson et a/., 2008, Gerschenfeld, 
2019 and Papukchiev et al., 2019). 


The simplest form of multiscale modelling is one-way coupling with a system code providing in- 
puts for a CFD model or a CFD model providing a correlation, response surface or Reduced Or- 
der Model (ROM) for a system code. However, this is often not appropriate for reactor scenarios 
and so two-way coupling is required. The process of code coupling can be complicated and can 
significantly impact the accuracy of a coupled simulation. Various coupling methods have been 
developed, which are discussed and classified by Zhang (2020), as follows: 


* Operation mode i.e. system and CFD codes run in serial or parallel. Parallel operation is 
more efficient and scalable, but is more complex to develop and maintain. 


* Coupling architecture describes how the codes are coupled (internal or external). Internal 
coupling means that one code is compiled inside the other i.e. single executable, while ex- 
ternal coupling means that the two codes are separate and communicate via input/output 
files, server-client (supervisor) or server-less (direct). Overall, the external server methods 
are considered the most efficient, flexible and maintainable. 


« Domain coupling is one of the critical parts of the coupling process, as it specifies which 
code is in charge of which computational domains. In domain overlapping, the system code 
simulates the whole domain and overlaps the CFD model, while domain composition has 
separate models with data transfer at the interfaces (Figure 3.1). Domain overlapping is 
more robust and efficient and means that the stand-alone and coupled system code models 
are identical, but requires more intensive code modifications. 


System Model 1 


Domain-Overlapping 


CFD Model 
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Figure 3.1: Domain overlapping and domain decomposition approaches (Zhang, 2020). 
This image is reproduced from an open access article under the terms of the 
Creative Commons Attribution License. 


¢ Field mapping refers to the method used to transfer data between the different spatial reso- 


33 of 84 


3.1.2 


Methodology 


lutions of the codes (CFD mesh is much finer than system code). This can be achieved using 
manual definition of mesh matching relations, user-developed toolkit for mesh matching or 
a third-party mesh-processing toolkit. A third-party toolkit is likely to be more accurate and 
flexible, but still requires significant code development to use. 


¢ Temporal coupling is needed to ensure that the multiscale code coupling is synchronous in 
time. This can be implemented using explicit, semi-explicit or implicit coupling. Convergence 
is generally better using a semi-implicit or implicit method, which should improve the accuracy 
and robustness. 


In general, it is important to ensure that the approaches are as consistent as possible between 
the two codes, especially if two-way coupling is being implemented, in areas such as material 
properties and boundary conditions. 


Multiphysics 


Multiphysics simulations involve coupling together multiple physical models or multiple simultane- 
ous physical phenomena. From a NTH perspective, this is predominantly associated with coupling 
thermal hydraulics analysis with neutronics and structural analysis (CSNI, 2015a). Since this set of 
technical volumes is primarily focused on thermal hydraulics analysis, only a brief introduction to 
multiphysics simulations is provided. 


Neutronics: The fuel cladding temperature and overall core behaviour depends on the balance 
between thermal hydraulics and neutronics. The thermal hydraulic behaviour depends on the heat 
release distribution due to radiation activity, while the neutronic behaviour is influenced by the 
fluid temperature, density and species concentration (e.g. boron). Most CFD software does not 
include any neutronics capability, while system/subchannel codes generally only include simple 
point kinetic neutronics models. Therefore, the coupling of CFD/system codes with more advanced 
(i.e. deterministic or stochastic transport) neutronics models will enable more accurate predictions 
of the fuel behaviour and fuel cladding temperature. An introduction to neutronics, and a discussion 
of thermal hydraulic/neutronic coupling, as well as guidance on more detailled references can be 
found in Voume 6. 


Fluid-Structure Interaction (FSI): Although temperature predictions are a key input to struc- 
tural integrity assessments (Section 2.4), coupling of CFD and structural codes is predominantly 
focused on Fluid-Structure Interaction (FSI), such as Flow Induced Vibration (FIV). This includes 
phenomena such as vortex-induced vibration, galloping/flutter, whirling, turbulence-induced vibra- 
tion and acoustic vibration (Kaneko et a/., 2014). FIV is a key issue for NPPs because primary 
circuit flows can induce vibrations that cause damage or wear to fuel assemblies, heat exchangers 
or pipework. Depending on the feedback between the structural displacement and flow field, this 
can be simulated using one-way or two-way coupling. Coupling with structural codes is available in 
most CFD codes, but care needs to be taken to resolve the vibration response, particularly at high 
frequencies, and this is the focus of the current CSNI CFD benchmark at the time of writing. 
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System and Subchannel Analysis 


System and subchannel codes are introduced in Volume 1 (Section 4.4.1) and their capability and 
limitations for modelling convection is discussed in Volume 3 (Section 3.1). 


System codes are a key design tool, allowing the performance of a whole system (or a group of 
systems) to be assessed in a wide range of scenarios much more quickly than would be possible 
with other analysis methods. By contrast, subchannel codes are used to perform reduced-order 
spatial analysis of complex components within a system (particularly the reactor core) for which 
the use of more detailed CFD would be too expensive or cumbersome. 


System and subchannel codes have the capability to model transient heat conduction within solid 
structures such as pipes and channel walls, albeit in a simplified manner compared to modelling the 
fluid flow. For many steady-state scenarios, or transient scenarios with largely isothermal mixing, 
these are generally reliable and robust. However, there are a number of scenarios involving CHT 
that are challenging to model using system and subchannel codes, some of which may present 
significant structural integrity issues. These physical scenarios include, but are not limited to: 


* Reflooding the core following a (large-break) LOCA where the overheated fuel rods are 
quenched to avoid fuel damage. 


« PTS may occur during overcooling events (e.g. ECCS water injection into cold legs and 
excessive heat removal by the secondary side). Natural circulation/flow stagnation is partic- 
ularly important for PTS as cold water is delivered directly to the downcomer, if the loop flow 
stagnates. 


* Temperature fluctuation may be observed due to turbulent mixing of two flow streams of 
different temperatures (e.g. in T-junctions). This phenomenon is often referred to as thermal 
striping, which can cause high cycle thermal stress and cracking. 


* Thermal stratification is a phenomenon where two flow streams of different temperature 
come into contact. The temperature difference causes colder (therefore more dense) water to 
occupy the bottom portion of a pipe, while warmer water flows over the colder water. Thermal 
stratification typically occurs in pipes where low flow velocity is expected. The circumferential 
temperature gradient can cause thermal stress which can lead to fatigue damage (Kim et al., 
1993). Thermal stratification can occur in pressuriser surge lines, feed water lines and spray 
lines, particularly during start-up and shut-down processes. 


« Flow reversal or counter-current flow can also occur in reactor components, such as steam 
generator tubes, upper plenum and core. This can lead to circumferential temperature gradi- 
ents that can cause thermal stress and fatigue damage. 


While system and subchannel codes may struggle to predict these scenarios, there are some 
special models implemented in these codes that can help improve the prediction of some of these 
phenomena. The following sections introduce system and subchannel codes, their application to 
CHT and highlight some of these particular aspects, but this should only be considered as a high 
level introduction; codes are typically accompanied by comprehensive user guides and validation 
information that should be consulted. 
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System Analysis 


In system codes, solid components can be modelled using heat structures (or heat slabs) that are 
represented using simple geometry (e.g. rectangle, cylinder and sphere). Heat structures allow the 
calculation of heat conduction in solids and can model: 


* Fuel pins/plates with nuclear or electrical heating. 
* Heat exchangers. 
* Heat transfer from pipe/vessel walls (e.g. environmental loss). 


Heat structures may be classified into two groups: active structures (e.g. fuel rods and heat ex- 
changers) and passive structures (e.g. pipe/vessel walls). While passive structures may not play 
an important role in a steady-state calculation, they may become significant heat sinks (or sources) 
during a transient scenario. Therefore, it is important to account for them in a simulation. Environ- 
mental heat loss may not be significant for full-scale NPPs and so may not be considered in a 
system code analysis. However, when modelling reduced-scale Integral Effect Test (IET) facilities, 
this heat loss can become significant, as the ratio of heat loss to core power becomes higher than 
that of a full-scale NPP. 


In general, both one-dimensional (i.e. radial) and two-dimensional (i.e. radial and axial) forms of 
heat conduction are available in system codes. The two-dimensional form is typically used for fuel 
components when reflood is expected. The predicted temperatures of solid components should be 
used with caution depending on the complexity of the component and the scenario being analysed 
in the system code. 


For components with a complex geometry and scenarios where the flow, pressure and coolant 
temperature conditions local to the component are changing rapidly, system codes may not accu- 
rately predict the true maximum or minimum temperature of the component. For components with 
a simple geometry in a scenario where the conditions local to the component are not changing 
rapidly (e.g. during steady-state operation or transient scenarios with largely isothermal mixing) 
system codes are more able to accurately predict solid temperatures. 


Nodalisation 


Users must normally specify the number of nodes in the radial direction of heat structures. As the 
computational cost can increase with increasing number of nodes, there is a trade-off between 
the solution accuracy and the number of nodes. In practice, a nodalisation sensitivity study may 
be necessary. The number of radial nodes typically varies between 2 and 20. If a solid structure 
(e.g. nuclear fuel rods) consists of different materials, a node should be placed at the interface 
between two materials. The mesh strategy can depend on the Biot number, Bi (Section 2.1.1.2). 
For example, the following recommendations are given for RELAP5-3D (Aumiller and Bettis, 2000): 


¢ For higher Bi (i.e. > 15), a graded nodalisation with smaller node intervals near the surface 
may be appropriate. 


* The recommended maximum surface node size can be expressed in terms of Bi as: 


Bmax = 1/(0.338Bi + 5.2) 
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where 6 = 0.56/L, and 6, is the distance from the surface node to the first internal solid 
node and L is the wall thickness. 


It should be noted that too small a nodalisation (< 10~* m) may lead to numerical oscillations 
in the solution. 


In the case of modelling reflood processes, the two-dimensional form of heat conduction should be 


considered. However, as the axial nodalisation is taken from the adjacent hydrodynamic volume, 


this can be too coarse to resolve the axial temperature gradient in the solid. In most system codes, 


an axial re-noding scheme is available to impose a finer nodalisation according to a set of criteria 


based on the local flow conditions. 


3.2.1.2 Models 


There are several special models available in most system codes, which are designed to deal with 


the scenarios described in Section 3.2 that can improve temperature predictions. These include: 


A reflood model that allows axial conduction, because this can play a significant role during 
the reflood process (Osakabe and Sudo, 1983). In addition, axial re-zoning or re-noding 
schemes may be available to resolve large axial gradients in temperature and/or heat flux. 


A Counter-Current Flow Limitation (CCFL) model may be considered when simulating re- 
flood and/or PTS at tie plates, steam generator U-tubes or hot legs. Without a CCFL model, 
the code may incorrectly predict the thermal hydraulic behaviour (e.g. zero liquid delivery). 
However, it should be noted that counter-current flow of a single phase in a horizontal pipe 
(e.g. pressuriser surge line) cannot be modelled in system codes due to its inherent vol- 
ume averaging technique. This therefore prevents the application of system codes to model 
thermal striping and thermal stratification in horizontal pipes. 


Some system codes may have a special model for thermal stratification in vertically oriented 
components (e.g. passive cooling pools) to limit the effect of numerical diffusion and improve 
the axial temperature profile prediction. 


Component specific models for nuclear fuel rods such as fuel deformation, gap conductance 
or metal-water reactions are available, however these models will contain approximations 
and so the user should review the impact of any simplification on the analysis. 


3.2.1.3 Limitations 


Other limitations of system codes to consider may include: 


System codes rely on a number of empirical correlations and a simplified set of governing 
equations which will have an impact on the accuracy of the temperature predictions, espe- 
cially during scenarios where the plant conditions are changing rapidly. 


System codes do not necessarily reach a nodalisation independent solution (by increasing 
the spatial resolution). 


An explicit representation of some details of a complex reactor geometry and structures may 
not be possible. 


Due to the simplified treatment of turbulence in system codes, prediction of turbulent thermal 
mixing is limited in system codes which will impact heat transfer predictions. 
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« While system codes have multi-dimensional flow modelling capability for components such 
as the reactor vessel, these are typically limited in their fidelity and accuracy. 


« Some system codes may have a tendency to predict lower Peak Cladding Temperatures 
(PCTs) and earlier quenching compared with experimental data (Choi and No, 2010, Li et a/., 
2016). This is often attributed to code deficiencies in the film boiling and interfacial drag 
models. 


These limitations should be considered carefully when undertaking system code analysis, and 
more detailed modelling should be considered if appropriate. For example, system codes could be 
used to provide initial and boundary conditions for higher fidelity analysis such as CFD or coupled 
with CFD codes to resolve complex local flow features (Section 3.1). 


Subchannel Analysis 


In subchannel codes, solid conductor models are available. In general, there are two types of mod- 
els: rods and generic conductors (e.g. vessel internal and external structures). Similar to system 
codes, solid conductors are represented as cylinders (solid or hollow) or plates (Moorthi et ai, 
2018). These allow subchannel codes to model vessel structures including: 


* Nuclear or electrically heated fuel pins. 
* Control rod guide tubes. 
« Support columns. 


In subchannel codes, an option to choose between one-dimensional (radial) and two-dimensional 
(radial and axial) forms of the heat conduction equation is available. Some codes such as CTF 
(Avramova and Cuervo, 2013) allow a three-dimensional (i.e. radial, axial and azimuthal) form 
of heat conduction. Multi-dimensional forms of heat conduction are recommended for cases with 
large temperature gradients in axial and/or azimuthal directions (e.g. reflood and quenching). 


Mesh Generation 


For cases with gradual heating/cooling, heat conduction in only the radial direction may be suffi- 
cient and the user needs to specify the number of radial points. While the number of radial mesh 
points must be sufficiently large to resolve the radial temperature profile, it can also increase the 
computational cost and so sensitivity studies may be required. 


In the case of nuclear fuel rods, a typical radial mesh size used in the fuel region is approximately 
10-3 m (Avramova and Cuervo, 2013). The multi-dimensional models may improve the code ac- 
curacy for scenarios involving large temperature gradients in the axial and/or azimuthal directions 
(e.g. during reflood). When axial conduction is selected, the axial mesh is typically taken from the 
adjacent fluid channel. 


The assumption of negligible azimuthal temperature gradients may be sufficient for many analyses, 
but it can be important in some cases. The azimuthal mesh depends on how a rod is connected 
to the subchannel region and tends to be divided into four azimuthal sections, following the chan- 
nel centred approach. Similar to system codes, subchannel codes can provide axial re-noding 
schemes for reflood calculation to impose a refined mesh according to a set of criteria based on 
local flow conditions. 
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There are several special models available in most subchannel codes, which are designed to deal 
with the scenarios described in Section 3.2 that will improve the temperature predictions. These 
include: 


The presence of spacer grids can influence flow conditions, such as cross-flow and flow 
blockage, and may improve heat transfer by introducing turbulence (Sugimoto and Muraov, 
1984). In most subchannel analyses, the effects of spacer grids are modelled using an effec- 
tive loss coefficient, although some subchannel codes also have additional models specif- 
ically for spacer grids including heat transfer enhancement. However, as the geometry of 
spacer grids is simplified in subchannel codes, the results of CFD simulations could be used 
to better capture the effect of spacer grids in subchannel codes (Salko et al., 2019). 


A turbulent mixing model is often available in subchannel codes to account for lateral turbu- 
lent exchange (Moorthi et a/., 2018). While this model is not specific to any phenomenon or 
scenario, it can impact the hydrodynamic volume solution which in turn can affect the solution 
of heat conductors. Users may be required to provide a mixing coefficient for single-phase 
and/or two-phase conditions and model sensitivity analysis may be required. However, the 
turbulent mixing model in most subchannel codes is simplified and applies the same turbulent 
diffusion coefficient to all mass, momentum and energy exchanges. 


As for system codes, special models for nuclear fuel rods (fuel deformation, gap conductance 
or metal-water reactions) may also be available in subchannel codes, and the impact of any 
simplification on the analysis should likewise be assessed. 


3.2.2.3 Limitations 


Other limitations of subchannel codes to consider may include: 


Subchannel codes depend on a number of empirical correlations and a simplified set of 
governing equations, so their suitability and validity needs to be understood for a particular 
application and software-specific guidance should be consulted. 


Subchannel codes are designed to model vertically-oriented components (e.g. reactor pres- 
sure vessel and core). Phenomena associated with horizontally-oriented components (e.g. 
thermal striping and thermal stratification) cannot be predicted. 


Only simplified heat conductor geometry is available. 


The fidelity of subchannel codes is lower than that of CFD and so may not be able to predict 
detailed local phenomena in channels. 


Despite these limitations, subchannel codes are useful as they allow for multi-dimensional mod- 
elling with reasonable computational time. The accuracy of subchannel codes may be improved by 
using CFD results to provide additional information (or coupling with CFD) to account for complex 
geometry effects. 
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CHT Approaches 


Common approaches to CHT analysis include uncoupled and coupled methods using a resolved 
solid, shell conduction or porous media zones. These are illustrated in Figure 3.2 and discussed 
in more detail in the following sections. A key reason for performing CHT analysis is to support 
structural integrity assessments (Section 2.4). Therefore, CHT analysis is often used to generate 
inputs for FEA models, enabling solid deformations, internal stresses and thermal fatigue to be 
predicted by structural engineers using tools that are familiar to them and well suited to these 
purposes. 


CFD Model 


CFD Model 


PACERS 


| Shell 
Structural Model CFD Model Pseudo Mesh CFD Model 
(a) Uncoupled (b) Coupled (c) Coupled (d) Shell (e) Porous 
Partitioned Monolithic Conduction Zone 


Figure 3.2: The main conjugate heat transfer modelling approaches. 


Uncoupled 


In this approach, the flow and fluid temperatures are predicted in a CFD model and the tempera- 
tures within solids are predicted in a separate thermal-structural FEA model. Information is passed 
one-way (from the CFD model to the FEA model) using data extracted from the flow solution. 
From a CFD perspective, the flow is predicted in the normal manner, but the following aspects are 
particularly significant: 


* Setting appropriate boundary conditions for the CFD model, given that the solids are not 
modelled (discussed further in Section 3.4.3). 


¢ Appropriate data will need to be extracted from the CFD solution in a form that can be used in 
the FEA model, often this is predicted HTCs and fluid temperatures or heat fluxes (discussed 
in Section 3.4.6). 


In view of the one-way information exchange (from fluid to solid), this method is generally only 
appropriate for situations where the thermal environment is changing so slowly that the fluids and 
solids can be considered to be in a quasi-steady-state, or the heat transfer coefficient is effectively 


40 of 84 


3.3.2 


Methodology 


independent of the solid temperatures (which may occur in forced convection flow). 


An uncoupled approach may be the simplest and computationally cheapest method for calculating 
the temperatures and stresses within components over long fault scenarios (hours or days). For 
example, the temperature variation in a weld in a pipe network or tube plate during start-up, shut- 
down or fault scenario. In this example, the inlet flow rate and temperature (calculated using a 
system level method) varies with time during the fault. The transient temperatures within the solid 
may be calculated by: 


¢ Calculating initial steady-state flow and temperatures within the component e.g. normal op- 
eration. 


* Identifying and calculating a number of quasi-steady-state flow conditions over the transient 
(e.g. when the inlet flow rate and/or fluid temperature changes). 


* For each condition, the fluid temperature and HTC is extracted using one of the methods 
described in Section 3.4.6. This could be a single value for a pipe or spatially varying value 
output at each face of the CFD model. 


* The fluid temperatures and HTCs are then imported into the thermal-structural model as a 
time-varying boundary condition, so that the solid temperatures can be calculated over the 
whole fault scenario. 


Coupled using Resolved Solid 


In a fully coupled CHT approach, the energy equation for temperature needs to be solved in both 
the fluid and solid domains. The solid temperatures from the CFD solution can then be directly 
mapped onto an FEA model to calculate the thermal stresses. The physical processes and solu- 
tions of the governing equations can be either solved in a single model (monolithic approach) or in 
separate fluid and solid models (partitioned approach). 


Monolithic or direct approach: The governing equations for both domains are conjugated by si- 

multaneously solving a single large set of equations (Patankar, 1980). A single generalised 
conservation expression is implemented for the fluid and solid domains. This ensures spatial 
and temporal energy conservation across the fluid and solid regions, and is simple to use 
within a single CFD code. However, if the temperatures are required as part of a structural 
analysis, then the spatial and temporal thermal field will need to be extracted at appropriate 
time steps and imported into the structural model. 
Although it is possible to separate the fluid and solid time steps, this approach often relies on 
explicit time integration for simplicity and is therefore penalised by small time stepping due to 
the most stringent phenomenon (Koren et al., 2017). It is also beneficial to consider whether 
a partitioned mesh for parallelised solution is balanced appropriately, because the solid cells 
in a partition will have significantly lower computational requirements than the fluid ones. 


Partitioned or iterative approach: This approach, described by Giles (1997), is based on a sep- 
arate solution of each set of equations for the fluid and solid, respectively. Each solution 
of the solid or fluid domains produces the necessary outputs along the fluid-solid interface 
which is then used as the boundary condition for the other. This iterative approach is widely 
used in finite difference, element and volume methods. By splitting the problem, an addi- 
tional common interface is created between the physical domains, through which energy is 
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exchanged, to allow the sub-domains to interact and yield a coupled solution. The use of a 
dedicated solver for the fluid and solid means that each can be optimised in terms of time 
step, algorithms or data structure. 


Most procedures in the literature impose a temperature boundary condition (Dirichlet type) 
on the fluid side of the interface and a pure (Neumann) or mixed (Robin) heat flux condition 
on the solid side. The Dirichlet-Robin coupling schemes are shown to offer many attractive 
features over the Dirichlet-Neumann schemes (El Khoury et a/., 2017, Errera et al., 2019), but 
requires selection of a suitable relaxation parameter which satisfies the stability and the fast 
convergence rate of the coupled problem (El Khoury et a/., 2017). The coupling is carried out 
so that both solvers compute temporal iterations separately until a given physical time that 
corresponds to the prescribed coupling time step (or ‘coupling period’). At these instants, 
data from both solvers are exchanged with each other in order to update their boundary 
conditions. With synchronized solvers, the fluid and solid solver time steps must fulfil the 
following constraint: 
Atep: = NsAts = Ne Ate 


where N, and N; are the number of iterations before exchanging data at the shared interface. 
Although partitioned approaches are flexible and enable direct coupling between the fluid and 
thermal-structural models, it is important to ensure that the interface is set up appropriately 
and that energy is conserved. This increases the user involvement and requires efficient 
dynamic coupling between the models. 

Koren et al. (2017) note that the burden of solving heat transfer in the solid after each LES 
time step is impractical in large-scale simulations with an order of magnitude increase in 
computational time compared to every 50 iterations giving a moderate 18% additional cost for 
the coupled simulation. However, since the varying wall heat flux and temperature depends 
on the size of the coupling time step, a large value might impact the accuracy of a LES 
prediction as the unsteady phenomena may be under-resolved at the fluid-solid interface. 
Therefore, it is important to check the impact of the coupling time step on the LES results. 


In both of these CHT approaches, the solid domain needs to be properly modelled and meshed. 
The parts of the computational domain containing fluids and solids are therefore both included 
in the computational mesh (Section 3.4.2). This can significantly add to the complexity and time 
required to build the model, especially if conformal meshes are used at the fluid-solid interface. 


A fully coupled approach is suitable for any combination of fluid and solid regions, including com- 
plex conduction paths and inter-connected fluid-solid regions. It also enables time-varying temper- 
atures to be calculated when the fluid and solid temperatures are strongly dependent on each other 
(e.g. mixed and natural convection, complex conduction paths and large Biot number). The main 
limitation of this approach is the need to properly resolve the timescales associated with the fluid 
phenomena, which often requires a small time step within the fluid domain. This means that long 
fault scenarios (Spanning hours or days) are often unachievable within realistic analysis timescales 
even for Reynolds-Averaged Navier-Stokes (RANS) modelling methods. The reduced timescales 
and larger meshes associated with LES models currently make this approach impractical for all but 
the simplest cases. 
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Coupled using Shell Conduction 


Creating full meshes for solids in order to resolve the temperature gradients within them can be 
onerous, yet for many geometries (such as pipes or ducting) the conduction behaviour is often 
quite straightforward. Therefore, some CFD software enables thin-wall or shell conduction to be 
modelled in the solver, by predicting the temperature variation within solids using a virtual (pseudo- 
mesh) that is generated by the CFD solver with a prescribed wall thickness. Although the imple- 
mentation and terminology will likely be code dependent, these terms are commonly described 
as: 


* Thin-wall conduction only accounts for heat transfer normal to the wall. 
¢ Shell conduction also includes planar conduction through the (virtual) solid. 


A wall/solid is typically represented by one or more layers which are added virtually by the compu- 
tational solver. Thus, there is no need to create the wall geometry or generate a solid mesh. The 
thermal resistance of the wall is accounted for by assigning appropriate material properties. This 
approach enables a number of solid layers to be accounted for at the same time (e.g. metal pipe 
and surrounding insulation). Within the same model, complex three-dimensional areas of geometry 
might be resolved using coupled methods with a mesh, and more straightforward areas modelled 
using shell conduction. 


External walls (fluid on one side): Shell conduction is ideally suited to applications with simple, 
uniform solid geometry (e.g. pipe walls, thin plates). This enables the thermal mass and conduction 
through the solid to be incorporated into the CFD solution. The main limitation is that the spatial 
temperature variation within the solid is not recorded as part of the solution, so this method is not 
appropriate if the solid temperatures are required for input into a thermal-structural calculation. 


Internal walls (fluid on both sides): Shell conduction can be used to simplify the model geom- 
etry by eliminating the need to mesh the solid and allow prism layers to be generated in the fluid on 
both sides more easily. As with external walls, this means that the thermal mass and conduction 
through the solid is accounted for, although it will also increase the flow area and volume of the 
fluid region to some extent. Therefore, the reduced mesh complexity needs to be balanced against 
the impact on flow accuracy depending on the purpose and requirements of the CFD analysis. For 
thin walls surrounded by large fluid domains, this may well be acceptable. 


Coupled using Porous Zones 


Most commercial CFD solvers offer the ability to specify porous zones. These can be applied to 
a wide variety of problems, including perforated plates, tube banks and pebble beds. This applies 
an additional pressure loss to the flow that is usually based on a user specified anisotropic loss 
coefficient. In laminar flows, the pressure drop is typically proportional to velocity (Darcy’s law), 
while in turbulent flows it is proportional to dynamic head (kinetic energy per unit volume of the 
flow). 


In most cases, the superficial velocity, based on the volumetric flow rate, is solved within the CFD 
software to ensure continuity of momentum and mass flow. Therefore, care is needed to ensure 
that the loss coefficient is appropriately derived, and the treatment of turbulence within the porous 
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zone is understood (by default, the turbulence generation and dissipation rate is often assumed to 
be unaffected by the solid component). The thermal equations can be solved using an equilibrium 
or non-equilibrium model: 


¢ Equilibrium model: The solid and fluid are locally assumed to be at the same temperature and 
are represented using effective material properties combining solid and fluid contributions, 
based on the porosity of the medium (percentage of fluid by volume). The thermal mass and 
conductivity of the solid components within the porous zone are taken into account using the 
zone porosity (Section 3.4.3.3). This enables an effective thermal conductivity of the porous 
zone to be calculated, and the thermal mass of the solid to be represented in the model. 


Non-equilibrium model: Separate temperatures are modelled for the solid and fluid zones 
with heat transfer between them. This approach is more appropriate for transient calculations 
because it enables changes in the solid temperature to ‘lag’ changes in the fluid tempera- 
ture, because of its representation of thermal mass and the finite rate of heat transfer. In 
addition, the effect of porosity on the transport velocity for energy and species in the flow can 
be included, even if the superficial velocity is solved for (although this might vary between 
codes). 


Porous zones provide a useful method of achieving the appropriate pressure drop and tempera- 
ture change in regions that are too complex to mesh explicitly. However, care must be taken to 
ensure that the impact associated with this approximation on outputs required from the analysis 
is understood. In some cases, a detailed CFD model of a complex component (e.g. tube bundle) 
could be generated and solved at a number of representative conditions. This could be used to 
develop pressure drop, turbulence and thermal correlations for an equivalent porous zone model. 
This would provide more confidence in the porous zone model for this application and enable the 
impact of the approximation to be quantified. 


CFD Analysis 


This section considers CHT analysis using CFD. It builds on the more general discussion in Vol- 
ume 1 (Section 4.5.3) and complements a number of more general topics given in Volume 3 (Sec- 
tion 3.2.2). General CHT approaches are considered in Section 3.3. As such, only areas particu- 
larly significant to modelling CHT using CFD are considered, such as approach, mesh generation, 
aspects of case definition and coupling. 


Guidance on industrial CFD computations (e.g. industrial guidelines such as ERCOFTAC, 2000, 
NAFEMS, 2019, and CSNI, 2015b, and books such as Versteeg and Malalasekera, 2007) are good 
general references, but their focus is not on providing detailed guidance for CHT analysis. 


CFD Approaches 


A number of CFD approaches are available, primarily resulting from different ways of predicting 
turbulence, and these are introduced in Volume 1 (Section 4.5.3). The application of these ap- 
proaches to CHT is introduced briefly below and discussed in more detail in the following sections. 
There are a range of approaches and numerous variations within each; for NTH, engineers must 
make informed judgements about the complexity and fidelity needed for their modelling work within 
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the context of a graded approach (Volume 1, Section 2.2.5) informed by their safety and economic 
case. More general NTH guidance is provided in D’Auria (2017) and CSNI (2015b). 


DNS: Direct Numerical Simulation (DNS) can resolve all length and time scales within a flow field. 
As discussed in Volume 1, it is considered highly accurate, but is very computationally expensive, 
so is predominantly used to study low Re or Gr flows in simple configurations. For CHT, DNS has 
been used to study basic scenarios such as fully-developed channel flows with steel walls, back- 
ward facing steps and T-junctions (Tiselj et a/., 2020) and jets in cross-flow (Wu ef a/., 2017). The 
extremely detailed data that DNS can provide has been used to develop fundamental understand- 
ing and refine the models used in other CFD approaches. However, in view of its limited industrial 
use, DNS is not considered further. 


LES: In LES, the Navier-Stokes equations are ‘filtered’, so that a proportion of the turbulence 
is resolved using the mesh and the remainder is modelled using a Sub-Grid-Scale (SGS) model. 
As discussed in Volume 1, LES is much less computationally expensive than DNS and potentially 
more accurate than RANS, but the computational costs (which scale with Re and Gr) are still high. 
Within the context of a graded approach for NTH, LES is most likely to be used in areas of high 
safety significance or commercial impact. For CHT, LES may be able to predict the fluctuating metal 
temperatures resulting from the turbulence within the flow field that can cause thermal fatigue (i.e. 
unsteady temperature fluctuations under steady flow conditions). Therefore LES has been used for 
assessments of primary circuit components, such as T-junctions (CSNI, 2011), where hot and cold 
streams mix. 


RANS: In RANS, the instantaneous flow variables in the Navier-Stokes equations are decom- 
posed into mean and fluctuating parts, so that the mean flow is resolved and the effects of turbu- 
lence are represented by a model (unsteadiness in the mean flow can be captured using URANS). 
As discussed in Volume 1, RANS approaches are used for the vast majority of industrial CFD (in- 
cluding for CHT analysis) and are often used to support more detailed analysis methods like LES 
or hybrid approaches. Within the context of a graded approach for NTH, RANS is therefore likely to 
be used in most situations where CFD is needed, with the complexity of the modelling work reflect- 
ing the complexity of the case and the safety significance and commercial impact associated with 
the work. Despite their well-known shortcomings (which are discussed in Volume 3, Section 2.2.4), 
RANS models can provide useful CHT results, such as steady and transient temperature fields 
within solid components. For analysis of solid temperatures over long plant transients, URANS is 
likely to be the most detailed CFD method that can be used with practical computational expense. 


Hybrid Methods: Hybrid methods combine RANS with aspects of LES within the same model, 
to take advantages of both. As discussed in Volume 1, there are a number of significantly different 
approaches available (and variants within these approaches). For CHT, hybrid methods may well 
be useful for modelling flows that feature unsteady thermal mixing, unsteady free-shear or highly 
separated flows away from walls. Resolving turbulent length scales in some areas of the domain is 
likely to increase the computational expense compared to RANS, but reduce it compared to LES. 
Like LES, the use of hybrid approaches is increasing in industry, and their use within a graded 
approach for NTH sits between the points made above for RANS and LES. 
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Mesh Generation 


The majority of industrial CFD calculations use a computational mesh on which the flow is solved. 
The mesh is within the computational domain and is fitted to the geometry of interest. Its purpose 
is to enable gradients in the flow variables to be resolved by the solver, and its quality can affect 
the results of a CHT calculation significantly. 


Computational Domain and Geometry 


These key aspects of CFD analysis are considered in Volume 1 (Section 4.5.2). Additional aspects 
that may be significant for CHT analysis include: 


* The chosen approach to modelling the solids (Section 3.3) is likely to affect what is included 
in the computational domain. Using a shell conduction approach is likely to reduce the size 
of the computational domain compared to what might be needed for a fully resolved solid 
approach. 


Where solids are being meshed, separate fluid and solid volumes are needed. These vol- 
umes should be coincident (so there are no gaps between them) and their features and 
curvature well aligned (errors may be created if the geometry is offset too much). Solid vol- 
umes arising from CAD models are likely to need a similar level of geometry ‘clean-up’ as 
the fluid volumes. 


The approach to defining boundary conditions within solids at the edge of the domain should 
be considered in the same way as fluids. It may be necessary to extend the domain (and 
hence add more mesh) to ensure the thermal boundary conditions do not have a distorting 
effect on the results of the analysis. 


Meshing Background 


General background on mesh types and quality is provided in Volume 1 (Section 4.5.2). The follow- 
ing sections consider meshing solid and fluid volumes in more detail; different meshes are required 
for LES and RANS work, so these are discussed separately (DNS is not included as it is little used 
in industry). Software specific guidance should also be considered to ensure the mesh is consis- 
tent with the CFD solver being used. Meshing often involves compromises, and the impact of these 
on the key results or figures of merit can be investigated using validation data, sensitivity studies 
and uncertainty quantification (Volume 1, Section 4.3). 


As noted in Volume 1, when a flow solution is available, the predicted flow should be visualised 
on top of the fluid mesh and reviewed by an experienced CFD user who can consider whether the 
mesh is capable of appropriately resolving gradients that may (or perhaps should) exist in the flow. 
The predicted temperatures should be visualised on top of the solid mesh in the same way. 


Software-specific guidance should also be consulted for guidance on mesh connections for CHT 
analysis. In general, meshes are preferably conformal for fluid to fluid mesh interfaces (i.e. every 
surface face in one volume has an equivalent surface face in the adjacent volume that it con- 
nects to, Volume 1, Section 4.5.2). For fluid to solid interfaces, meshes do not generally need to 
be conformal. However, the temperatures and heat fluxes in both the fluids and solids should be 
visualised together around the interfaces, to check for gaps, offsets and interpolation errors (partic- 
ularly around any large differences in mesh sizing at an interface between fluid and solid volumes). 
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Additional checks such as reporting areas of interfaces and ‘unmapped nodes’ may be helpful to 
avoid errors. 


Meshing Solid Volumes 


Meshing solid volumes is often more straightforward than meshing fluid volumes (since gradients 
in flow fields do not need to be considered). However, it is still important that the mesh can capture 
the temperature gradients within the solid, and assessing Bi may help to improve understanding 
of these (Section 2.1.1). 


If large temperature gradients are expected in the solid or there are small geometrical features 
containing large temperature gradients the mesh should generally be sufficiently refined to capture 
these. For thin solid components (e.g. pipe walls) structured hexahedral or prism cells (rather than 
unstructured mesh) are recommended through the thickness to capture the thermal gradient. This 
may be particularly important if the temperatures of the solid are used to perform structural integrity 
assessments. If this is the case, it is good practice to discuss the work with the structural engineer 
at the outset, so that the mesh and data exchange process can be designed correctly from the 
start. 


For unsteady CHT analysis, a prismatic layer is recommended in the solid with a fine mesh at 
the fluid-solid interface to resolve the transient temperature profile in the solid. To demonstrate 
this, consider a constant surface heat flux into a solid at a uniform starting temperature; the initial 
temperature gradient at the wall boundary condition is fixed by the heat flux and the solid surface 
(wall) temperature, 7,,, depends on the size of the first cell (7, is the cell centre temperature of the 
first cell in the solid). 


k gAx 
= —(Ty-Ts), “. Tw = —— 
G=7,\ ) 5 


This means that if the first cell is too large, the wall surface temperature is over-predicted. In 


+ Ts 


practice, fixed heat flux conditions arising from an adjacent fluid do not occur, but the wall resolution 
does affect the predicted heat flux and interface temperature between the fluid and solid regions. 
Therefore, a thin first cell layer is needed in the solid to correctly predict the surface temperature, 
and since the thermal gradient in the solid and thermal response reduce with distance from the 
wall, a gradual increase in cell size is appropriate. 


The benefit of using a prismatic layer can be easily visualised by solving a simple plate conduction 
thermal transient (Figure 3.3). A thick (L = 2cm) stainless steel plate is initially at 600 K with a con- 
vection thermal boundary condition applied (HTC = 1000 W/m? K, T+ = 620K). The time-varying 
surface temperature and heat flux for an inflation layer mesh and simple coarse mesh (Figure 3.3c) 
show that the surface temperature is initially 3K too high for the coarse mesh, which reduces the 
heat flux into the plate. The inaccuracy of the surface temperature and through-thickness temper- 
ature profile (Figure 3.3d) for the coarse mesh reduces with time as the transient progresses. This 
highlights the importance of appropriately resolving the mesh in solids for unsteady CHT analy- 
sis or accurate surface temperature predictions. This is particularly relevant to applications that 
involve rapid fluctuations in the fluid temperature, such as thermal striping phenomena, PTS or 
ECCS injection. 
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Figure 3.3: Impact of mesh density on transient temperatures in solid. 


Meshing Fluid Volumes: LES 


In LES, the filtering (between what turbulence is resolved and what is modelled) is usually related 
to the mesh resolution. Meshes for LES must therefore be sufficiently fine to enable appropriate 
resolution of turbulence length scales in three dimensions. It is important to ensure that a wide 
enough range of turbulence length scales is resolved, otherwise the assumptions LES is based on 
may not be valid. This range is generally considered to be that at least 80% of the turbulent kinetic 
energy in the flow is resolved and the remainder is accounted for by the SGS model (Pope, 2000 
and Wagner et a/., 2007). 


LES calculations are inherently transient, so a time step is required by the solver. It is important 
to note that this time step must be appropriate for the flow field, and also consistent with the 
mesh size (Volume 3, Section 3.2.5). Higher Re or Gr flows often need a finer mesh and therefore 
generally use smaller time steps. Since the appropriate mesh size to use depends on the flow 
field, it is normal to start LES work by performing an initial RANS calculation. This gives an initial 
solution and develops understanding of the flow under study. It also enables the required LES 
mesh dimensions to be estimated based on the RANS turbulence model predictions by calculating 
the integral length scale and Taylor microscale values, together with a representative cell length. 
Software-specific guidance should be checked because definitions may vary (particularly for k-w 
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models), but equations normally have a form similar to: 


k3/2 ki/2 
Integral length scale (m): /g = —— = Gu’ where C,, = 0.09 (Menter, 2015) 
E ww 


|. k ; on ; 
Taylor microscale (m): > = 4/10v—, where v is the kinematic viscosity 
& 


Cell length (m): A = V/3 for modest aspect ratios, where V..) is the cell volume 


cell 


Typically, 10 to 20 mesh cells are needed across /p to resolve 80% of the turbulent kinetic energy, 
so the suggested cell size is typically taken to be A = max(A, /g/10) for full LES resolution (further 
detail is available in Addad et al., 2008). This can be checked by plotting /9/A < 10 and »/A < 1 
to identify areas that are under-resolved. Once the LES calculation is complete, its success in 
capturing the turbulence length scales can be checked by calculating the amount of turbulent kinetic 
energy that is resolved as a proportion of the total turbulent kinetic energy (Pope, 2000). The 
proportion of turbulent kinetic energy that is resolved in a LES solution can be estimated by: 


o = kres 
kres — kegs 
1 LL . 
Kres = - (Uzms + Vims + Wems) and kegs (4 ) 


Where k,., is the resolved turbulent kinetic energy, kg; is the time-averaged sub-grid scale (mod- 
elled) turbulent kinetic energy, Urms, Vrms ANd Wms are the Root Mean Square (RMS) fluctuations 
in the three components of velocity, jus; is the sub-grid scale turbulent viscosity, p is the density, 
C is the Smagorinsky coefficient and A is the sub-grid scale filter length. A value of o > 0.8 is 
considered to represent a good LES solution (Pope, 2000). It is worth noting that this definition 
of the proportion of turbulent kinetic energy resolved does not include the turbulent kinetic energy 
numerical discretisation error (Anum), which according to Celik et al. (2009) may be of similar mag- 
nitude to k.,;. Further guidance on the methods used to assess the quality and reliability of LES 
models are detailed in Meyers et al. (2008). 


LES is relatively intolerant to large cell aspect ratios and rapid changes in mesh size (because 
the turbulent viscosity is usually directly related to the local grid size). For boundary layer flows, a 
small near-wall mesh sizing (yt < 1) is required to model the small eddies in the inner-most part 
of the boundary layer (which are broadly comparable to wall distance in size). This therefore leads 
to similarly small mesh sizes in both streamwise and spanwise directions. An LES calculation that 
has properly resolved the near-wall scales is usually referred to as ‘wall-resolved’ LES. For most 
engineering flows, wall resolving leads to very onerous mesh requirements’. 


The near-wall region is often significant because the flow here directly affects heat transfer (Sec- 
tion 2.1.2). There have been attempts at developing wall functions for LES (see, for example, 
Piomelli, 2008), but these may be more restrictive than RANS counterparts since they need as- 
sumptions about how velocity fluctuations vary near the wall. However, in industrial cases it is 
common to use hybrid methods such as Wall Modeled Large Eddy Simulation (WMLES) to reduce 


For Reynolds number Re > 10°, this becomes prohibitive quite quickly; Piomelli (2014) estimates that over 90% of the grid 
points may be needed to resolve less than 10% of the computational domain. 
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the expense of modelling near-wall flows (or high Re flows may be unresolved near walls). 


It is noted that mesh sensitivity studies need care for LES, since finer meshes usually result in a 
greater proportion of the turbulent length scales being resolved (as opposed to being modelled) and 
the approach tends towards DNS. It is therefore not possible to achieve true ‘mesh independence’ 
in the flow fields in the same manner as is considered desirable for RANS approaches. Therefore, 
the quality of an LES simulation needs to be assessed based on cell aspect ratio and percentage of 
turbulent kinetic energy that is resolved, and potentially the convergence of time-averaged integral 
quantities (forces or heat fluxes on surfaces), or the frequency spectrum of fluctuating quantities of 
interest. 


Meshing Fluid Volumes: RANS 


For RANS models it is desirable to obtain mesh-independent solutions, so it is normal to use a 
mesh sensitivity study to understand the influence of the mesh on the solution. Ideally, at least three 
meshes of increasing resolution might be tested and compared using flow quantities of interest, 
ideally by doubling the mesh in all dimensions (increasing the mesh size by a factor of eight). 
However, due to computational expense, the impact of mesh resolution might be studied in key 
areas within the full computational domain, or a part of the computational domain might be meshed 
separately to study. This is discussed in more detail in Volume 4 (Confidence and Uncertainty). This 
section provides some guidance that may assist in developing a mesh-independent solution. 


Near Walls: As introduced in Volume 1, a key aspect of RANS methods is the approach to mod- 
elling flow near walls. Methods introduced in Volume 1 include wall resolving (low-Re, fine mesh, 
small yt), standard wall functions (high-Re, coarser mesh, larger y*), enhanced wall functions 
(CFD code tailors approach to local y*) and advanced wall functions (more complex models cater- 
ing for non-equilibrium conditions), although the terms used and detailed implementations may 
vary between CFD software. Each approach has different requirements for the local mesh spacing 
near the wall, which is normally defined in meshing software before a flow solution is available. 
Software specific guidance should be consulted, but enhanced wall functions are often used for 
industrial flows, so that: 


1. The CFD user decides where in the domain more detailed modelling of flows near walls is 
valuable (e.g. where there is significant surface heat transfer, or separating flows are ex- 
pected) and designs the mesh in these areas to resolve the viscous sublayers of turbulent 
boundary layers (which occur between 0 < y* <5, Schlichting and Gersten, 2017, Chapter 
17), consistent with a wall-resolving approach. To achieve this, the mesh is generally devel- 
oped so that the nearest point to the wall at which the solver will compute flow variables is at 
yt < 2, ideally yt = 1. The mesh should then expand smoothly away from the wall. 


2. In the other areas (where the CFD user considers detailed wall modelling less valuable, such 
as pipework or volumes where the details of the flow have less impact on the key results) 
the mesh is designed to only capture the fully turbulent region of boundary layers, consistent 
with a standard wall function approach. To achieve this, the mesh is generally developed so 
that the nearest point to the wall at which the solver will compute flow variables is at yt ~ 30 
to 50 (the CFD solver can then apply scaled wall functions in these areas), but higher values 
of 100 of more may be permissible, depending on the situation. 
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Whatever approach is taken, the mesh often needs to be adjusted when an initial flow solution is 
available. If the mesh is poorly designed however, these initial flow solutions will be misleading, 
and more iterations will be needed to obtain an appropriate mesh. The following guidance may 
help reduce the amount of iteration needed. It may also be helpful to develop meshes using small 
test cases on parts of the geometry before completing the full mesh. 


First Cell Height: To develop a near-wall mesh, the target y* value at a given point on a wall is 
turned into a distance between the wall and the nearest point to the wall at which the solver will 
compute flow variables (y). However, this can only be done with some understanding of the flow 
field, since the time-averaged wall shear stress (7) is needed: 


y= ee = ee 
Ur [ee 
p 


In practice, y is often set as a constant over a given surface (i.e. is not changed in response to the 
local flow field) and therefore: 


V L 1 
at aa it 
Y = Vearget = Yearget = Yearget 
a pce Ee fe OE ee) 
p PUN 2 Le 2 


Re/Lre is the Re per unit length, which can be assessed using estimated flow-field quantities 


(the lengths should cancel). pU can be calculated from a mass flow rate where appropriate (WV = 
pAcsU). C+ is the local skin (or ‘Fanning’ friction coefficient commonly used in external flow, which 
is related to the pipe (or ‘Darcy’, or ‘Moody’) friction coefficient (f) commonly used in internal flow 


as: 
f Tw 
4 3 pU? 


Ct or f can be estimated using a resistance diagram (or ‘Moody chart’ for pipes) or a correlation 
(Schlichting and Gersten, 2017; Rohsenow et al., 1998; Kakag et a/., 1987; Miller, 2009; Idelcik 
and Ginevskil, 2007). Example correlations for incompressible external flow over a smooth flat 
plate include: 


t = (2logy9(Rex) — 0.65)~*3 for turbulent flow (after Schlichting) 


Cl. = 0.664(Re,)~? for laminar flow 


Where Re, is based on the streamwise length x and free-stream velocity (laminar and turbulent 
flow are discussed in Volume 3, Section 2.2). Example correlations for fully developed flow through 
a straight rough pipe with circular cross-section include: 


—2 
5.74 

f =0.25 Ce (35 + as )| for turbulent flow (Swamee-Jain equation) 
. D 


f= spe. for laminar flow 
Rep 


Where Rep is based on the pipe diameter and superficial velocity (Section 2.1.2.3) and « is the 
surface roughness. The Swamee-Jain equation is one approximation of the implicit Colebrook- 
White equation (a number of other similar approximations are available). 
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As these diagrams and correlations are generally for simple geometry like flat plates or straight 
pipes, the predicted C; or f values will need to be increased to account for local acceleration of the 
flow over other geometries (which generally increase the wall shear stress). Engineering judgement 
is needed for this (as a rough starting point for external flow, a factor of 2 for streamlined shapes 
and 4 for shapes like cylinders may be useful). 


Once y has been predicted, the height of the first cell away from the wall (the ‘first cell height’) can 
be specified in the meshing software. For vertex centred numerical schemes the First Cell Height 
(FCH) is normally equal to y, while for cell centred schemes it is normally equal to 2y. 


Near Wall Prismatic Mesh: As noted in Section 2.1.2, external flows often feature boundary 
layers near walls surrounded by a free-stream flow, while internal flows may feature significant 
gradients in flow variables across the domain. In either case however, the near-wall regions are 
likely to feature large gradients in flow variables normal to the wall, and smaller gradients along the 
wall. These areas are often meshed using a prismatic mesh (Volume 1, Section 4.5.2). As well as 
the FCH, the thickness of the prismatic mesh and the number and distribution of mesh cells across 
this thickness should be appropriate for the flow (Figure 3.4). 


Unstructured Mesh 


Boundary Layer Thickness 


Prismatic Mesh Height 


t First Cell Height 


Figure 3.4: Illustration of FCH and near-wall prismatic mesh containing a boundary layer. 


In general, the prismatic mesh should be designed so it is thick enough to accommodate the full 
thickness of the boundary layer. Boundary layer thickness (6) for external flows can be predicted 
using a suitable correlation (Schlichting and Gersten, 2017). Example correlations for incompress- 
ible external flow over a flat plate include (Anderson, 2017): 


6= snl for turbulent flow 
Rel/® 
5.0x . 
6= Rel? for laminar flow 


Where Re, is evaluated at streamwise length x and free-stream velocity U.. As for the prediction 
of FCH, correlations are likely to be for simple geometry like flat plates, so predicted thicknesses 
will need to be increased to account for local deceleration of the flow (which generally increases 
boundary layer thickness). Similar factors to those used for estimating FCH above may provide 
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a rough starting point. For internal flows, the prismatic and adjoining mesh should be designed 
together to capture the expected flow gradients. 


An appropriate number of cells should be placed across the boundary layer thickness, and these 
should be distributed correctly. The distribution is normally set using an expansion factor from 
the FCH out into the flow. Since the height of the boundary layer is likely to be growing within 
the prismatic mesh along the flow direction, either a) the number of cells across the height of 
the prismatic mesh should be larger than the number of cells required across the thickness of 
the boundary later, or b) the height of the prismatic mesh should be tailored to the thickness of 
boundary layer as it grows (by tailoring the expansion factor). The former approach is normally 
used in industrial cases, with the latter approach normally only applied for applications like aircraft 
wing or turbomachinery aerodynamics. 


In practice, the prismatic mesh is normally defined using the FCH, a linear growth rate and the total 
number of cells. The prismatic mesh thickness is therefore not specified directly, and some iteration 
is needed to identify appropriate settings. Software specific guidance should be considered, but in 
general at least 10 to 20 mesh cells across the boundary layer thickness are often recommended, 
and 30 may be used in more detailed models (15 cells are often used as a minimum to reduce the 
risk of dropping below 10). However, this will depend on the wall treatment approach that is being 
used (i.e. wall-resolving or wall function method). 


The thickness and number of cells across the prismatic mesh can be checked (using an initial 
flow solution) by visualising the turbulent viscosity across the boundary layer. This quantity should 
have a maximum in the middle of the boundary layer, so the thickness of the boundary layer can 
be estimated by doubling the distance between the wall and where this maximum occurs, and the 
number of cells across this distance can also be estimated. However, if the boundary layer grows 
out of the prismatic mesh and into a coarser or more dissipative region of adjacent mesh, the 
flow field may not develop properly and boundary layer thicknesses may be effectively limited. If 
this occurs, it may be difficult to know how thick the boundary layer mesh should have been and 
hence what expansion ratios should have been used making it difficult to know what adjustments to 
make to the mesh to predict the flow correctly. Allowing a greater thickness of prismatic mesh and 
reducing it later if needed may be helpful. It may not be practical or possible to create a prismatic 
boundary layer mesh that is large enough in all situations; in these cases (and in general) the 
size and quality of cells (which may be unstructured) immediately outside of the prism layer that 
transition to the wider mesh require care and attention, particularly to avoid large jumps in cell size 
and wall-normal resolution. The field values for the solution in these areas should be scrutinised. 


Since temperature gradients are important for heat transfer, thermal boundary layers also need 
to be considered. This is important regardless of the wall boundary condition (adiabatic, specified 
temperature, specified heat flux, specified heat transfer coefficient etc). If Pr is close to unity, the 
Reynolds analogy can be used because the momentum and thermal boundary layers are likely 
to be similar (i.e. C;/2 = St). However, if Pr is not close to unity, the thermal and momentum 
boundary layers can be quite different; more detail is available in Volume 5. 
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Away From Walls: While the focus of this volume is on heat transfer to surfaces, it is also impor- 
tant that gradients in the wider flow field are appropriately resolved (especially for internal flows). 
As for boundary layers, enough mesh should be provided across areas of sheared flow (10 to 20 
mesh cells across a shear flow are often used, but software specific guidance should be checked). 
It is also important to control the aspect ratio of cells near shearing flows. This is particularly true 
for natural convection cases, so these aspects are discussed further in Volume 3 (Section 3.2.3). 


Case Definition 


Once a computational mesh has been developed, the CFD case can be configured by defining 
physical information and solver settings. Setting fluid boundary conditions, fluid material properties 
and spatial and temporal discretisation are considered in Volume 3 (Section 3.2.4) as these are 
especially relevant to natural convection. 


This section discusses additional aspects relevant to CHT analysis, considering boundary condi- 
tions, thermal radiation, the use of porous zones, and aspects associated with fluid to solid inter- 
faces. 


Boundary Conditions 


Boundary conditions define the problem being solved, so if they contain inaccurate information or 
are applied incorrectly, it is likely the results will also be inaccurate. The application of boundary 
conditions in CFD tools varies, so software specific documentation should be consulted. For CHT 
analysis, the below aspects are particularly significant: 


¢ Inflow and outflow boundary conditions should ideally be placed sufficiently far from areas 
of interest to avoid unwanted interactions. If an inflow is fully developed, suitable profiles for 
flow variables could be applied from literature sources, or created (e.g. by using a sub-model 
or extending the domain). The turbulence and temperature conditions for backflow at outlets 
may need to be carefully chosen, or backflow prevented by solver settings. Inlets and outlets 
also can have radiation settings applied so that the interior of the domain can ‘see’ a specified 
external temperature. 


¢ Wall boundary conditions are important because they affect the pressure drop and heat 
transfer in the fluid domain. For CHT analysis, the surface roughness and thickness of de- 
posits can have a significant impact on the shear stress, and hence heat transfer into the 
fluid (Section 2.2). The surface roughness should be estimated based on the best available 
knowledge, and potentially reviewed as part of a sensitivity assessment if through life con- 
siderations are important. Different CFD tools have different approaches and requirements 
for the relative size of the first cell height and the applied surface roughness; the treatment 
of roughness impacts the accuracy of the near-wall models and type of near-wall modelling 
approach that can be used. 


« External heat transfer is often a key requirement in CHT analysis and is applied to the ex- 
ternal solid surface of the model. This needs to be considered in the context of the analysis 
to assess whether the temperatures or thermal gradients are conservative, because there 
is often a lot of uncertainty and variation in the external conditions. The following methods 
are commonly used, and the values should be potentially reviewed as part of a sensitivity 
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assessment: 
— Adiabatic: A zero heat flux boundary condition is often used for surfaces that are cov- 
ered by thermal insulation or to maximise temperature predictions when the external 
heat loss is considered insignificant. 


— Constant heat flux: This applies a constant heat flux (W/m?) to the surface that is 
independent of the bulk temperature and flow conditions, such as electric heating ele- 
ments or heat released as a result of radioactive decay or the fission process. 


— Constant temperature: This fixes the temperature of a surface to a specified value and 
can be used to represent high external heat transfer conditions, such as pump driven 
water cooling, or boiling at a known saturation temperature. 


— Applied HTC and ambient temperature: This is the most flexible and generic ap- 
proach to applying external thermal boundary conditions, and involves specifying con- 
vection (HTC and bulk temperature) and thermal radiation (far-field temperature and 
emissivity) conditions. Therefore, the external heat transfer can vary and depends on 
the local wall temperature, which enables a more realistic representation of the exter- 
nal heat loss. For example, this can be used for forced or natural convection to air by 
using empirical correlations for HTC, such as a vertical or horizontal heated plate (In- 
cropera et a/., 2011), although the limitations and uncertainty associated with empirical 
correlations should be borne in mind (Section 2.1.2.3). 


¢ For thin-walls or shells (Section 3.3.3), the thickness and solid material properties need to 
be specified to enable the thermal gradient and thermal mass associated with the solid to be 
included in the model. This includes contact thermal resistance between two solid materials 
(Section 2.1.1.1), if this is considered significant. 


When porous zones are included (Section 3.3.4), the porosity (open volume fraction), pres- 
sure loss and effective material properties need to be defined (Section 3.4.3.3). In addition, 
the impact on turbulence generation and appropriate representation of the porous region to 
the physical geometry needs to be considered. 


For uncoupled or coupled partitioned models (Section 3.3), the interfaces between the CFD 
and FEA models need to be carefully defined and appropriate boundary conditions used to 
enforce continuity of temperature and heat flux across the two regions. This will depend on 
the purpose of the analysis and software-specific recommended approach. 


For coupled monolithic models (Section 3.3.2), the coupling between the fluid and solid re- 
gions is usually simple, as it is addressed within a single code. For non-conformal interfaces, 
it is important to ensure that heat transfer across the interface is enabled and large variations 
in the mesh resolution across the interface do not impact the accuracy of the temperature 
predictions. 
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Thermal Radiation 


Radiative heat transfer is complex (Section 2.1.3), and it is generally only applied in a simplified 
and approximate manner in NTH analyses. The need to model thermal radiation depends on the 
properties of the fluid: 


« Water and liquid metal: Radiative heat transfer is not important because thermal radiation 
does not penetrate any appreciable distance before it is absorbed. 


¢ Air and most simple gases (e.g. helium and argon): Radiative heat transfer occurs between 
surfaces and they are considered non-participating media (neither absorb, emit or scatter 
thermal radiation at the IR wavelengths of interest), and variations of surface properties with 
temperature and wavelength are usually ignored. 


¢ Water vapour, COz and soot: For large volumes and significant concentrations, they emit, 
absorb and scatter thermal radiation, and so may need to be considered as a participating 
medium. 


* Molten salts: These are semi-transparent media that absorb the different wavelengths emit- 
ted from a hot surface to different degrees according to its absorption spectrum and salt 
composition. This is discussed in detail in Volume 6. 


Thermal radiation becomes significant when the radiant heat flux, grag = €o(Tf 


max Tints is 
large compared to the heat flux due to conduction or convection, which typically occurs at high 


temperatures. 


However, radiative heat transfer can significantly impact the overall flow and thermal solution, and 
so should only be neglected after careful consideration. For example, Martyushev ef al. (2014) 
investigated the unsteady regimes of convective-radiative heat transfer in a cubic enclosure with 
finitely thick heat-conducting walls in the presence of a constant-temperature energy source. The 
dominant mechanism of transfer of energy in this system is natural convection, but neglecting 
the radiative component can lead to considerable variations in integral characteristics and in the 
structure of the velocity and temperature fields. 


There are a number of standard thermal radiation models available in most CFD software with 
different levels of complexity and medium participation. The most common models include: 


P,: The Radiative Transfer Equation (RTE) is simplified by the P, model to an approximate form 
that is fast to solve. It is mostly recommended for ‘optically thick’ (participating medium) ap- 
plications, and does not include significant detail on directionality, so is not recommended for 
complex geometry or transparent media. The P, model can include the effect of absorption 
and particulates within the medium with properties that can vary with radiation wavelength 
i.e. non-gray and gray radiation. It assumes that the reflection of incident radiation at a sur- 
face is isotropic (diffuse). 


Rosseland: When a medium is highly absorbing or scattering, then the Rosseland diffusion equa- 
tion can be used as an approximation for radiation transport. The diffusion approximation 
does not retain any directionality in radiation and heat transfer appears as an additional ra- 
diative conductivity (that depends on T°) added into the energy equation, so no additional 
radiation transport equations are solved. Because the diffusion approximation is not valid 
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near to a surface, ‘temperature jump’ (discontinuous) boundary conditions are needed when 
using this model. 


Surface-to-Surface (S2S): The S2S radiation model is useful for modelling radiative heat transfer 
within an enclosure with non-participating media. It calculates the view factor (Section 2.1.3.3) 
for each face or cluster of faces to determine the energy transfer between two surfaces based 
on their size, separation distance and orientation (and does not consider the medium be- 
tween them). The view factor calculation is usually undertaken as a pre-processing step, 
making the model faster to solve. The S2S model usually assumes gray radiation and that all 
surfaces are diffuse. 


Discrete Ordinates (DO): The DO model has the capability to solve radiative heat transfer for par- 
ticipating and non-participating media, including semi-transparent walls (Murthy and Mathur, 
1998). Therefore, it can include the effect of absorption, scattering and particulates within 
the medium for gray and non-gray radiation. It solves the RTE as a transport equation for 
radiation intensity for a finite number of discrete solid angles, and so the accuracy of the 
solution depends on the angular discretization. If the discretization is too coarse, it can lead 
to physically unrealistic ‘oumps’ and artfeacts in the incident radiation field. The main disad- 
vantage of the DO model is that it can be computationally expensive, and so increases the 
solution time. 


Monte Carlo (MC): The MC model simulates radiative heat transfer by tracking (‘ray-tracing’) a 
large number of representative photon bundles (histories) through the domain. It is a general- 
purpose method and also allows absorption, scattering and particulates within the medium 
for gray and non-gray radiation. MC simulations are also a balance between accuracy and 
computational costs, where improved accuracy is possible by increasing the number of pho- 
ton histories tracked, making the simulation computationally more expensive. Some MC 
implementations have an S2S option for non-participating media, where the photons are 
tracked from wall to wall without interacting with the medium they are passing through. 


For non-participating media, the key parameter is the emissivity of the surfaces in contact with the 
radiative heat transfer enclosure (Section 2.1.3.2). For participating media, it is also necessary to 
specify the absorption coefficient, scattering coefficient and refractive index for the fluid within the 
enclosure. 


Porous Zones 


As discussed in Section 3.3.4, porous zones can be used in CFD models to represent solid and 
fluid regions that are too complex to model explicitly, such as packed beds, perforated plates, filters 
and tube banks. Porous media models use an empirically calculated flow resistance that is applied 
to the porous zone. This is usually based on the superficial velocity (the velocity that the fluid would 
have through the zone without the solid component being present) and depends on the physical 
structure/arrangement of the porous medium and porosity (proportion of fluid by volume, ¢). 


The use of a superficial velocity enables the bulk pressure loss through the porous zone to be 
calculated, but it does not resolve the increased physical velocities that occur, which limits the ac- 
curacy of the model. More advanced porous media models based on the physical velocity do exist, 
but their assumptions and appropriateness for each application needs to be carefully considered. 
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The pressure drop in porous zones is usually represented by a viscous loss (laminar flow) and an 
inertial loss (turbulent flow). In laminar flow through porous media (e.g. packed beds), the pressure 
drop is typically proportional to the velocity and tends towards Darcy’s Law, which depends on the 
permeability of the medium, k (m?). 


b 


Laminar flow, AP = rade 


In turbulent flow through porous media (e.g. tube banks or perforated sheets), the pressure drop is 
typically proportional to the dynamic pressure and depends on the loss coefficient per metre, K/L 
(1/m). Typical values of loss coefficients for standard geometrical arrangements are provided in 
Miller (2009), and are discussed in Volume 3 (Section 2.2.2). In some cases, it may be appropriate 
to use anisotropic values of pressure loss based on the geometry (such as a rod bundle where 
losses in cross-flow are different to those in the axial direction). 


K1 
Turbulent flow, AP = The 


The expression for the momentum sink that is added to the Navier-Stokes equations in CFD codes 
that contains both laminar and turbulent parts is known as the Darcy-Forchheimer equation, or 
sometimes the Darcy equation with the Forchheimer (turbulent) correction. 


By default, most porous media models neglect the presence of the solid for turbulence generation 
or dissipation, which limits the accuracy of the model, although some modifications to the k-« 
equations have been proposed to account for the effects of porosity on turbulence (Al-Aabidy et ai., 
2020). For porous beds, turbulence will be reduced by the medium, while for cross-flow in tube 
banks turbulence is likely to be increased. Better understanding and accuracy of the porous media 
model could be achieved by using a CFD model to explicitly resolve the fluid and solid components 
in order to determine more accurate correlations for pressure drop and turbulence generation. 


Examples of the derivation and application of porous models are shown in Study B: Fuel Assembly 
CFD and UQ for a Molten Salt Reactor and Study C: Reactor Scale CFD for Decay Heat Removal 
in a Lead-cooled Fast Reactor. 


Effective Material Properties: For porous zones, it is necessary to calculate the effective ma- 
terial properties of the combined fluid and solid components. The effective properties are often 
automatically calculated in analysis software based on the porosity, but the method should be 
checked to confirm that it is appropriate for each application. 


The effective volumetric heat capacity (c,), or thermal mass, of a porous region is usually calcu- 
lated based on the overall porosity, as follows: 


p Pr 
Peff = (1 = ?)Ps =F pr and Cp,eff = (1 ¢) 2 Cp,s + ) Cp, 
Peff Peff 


J. (Cpe = (1 — $)(PCp)s + $(PCp)¢ 


However, thermal conductivity is harder to calculate as it depends on the structure/composition of 
the porous media. The effective thermal conductivity (kes) is defined as the thermal conductivity of 
the medium with no fluid motion. For general porous structures, the true value will lie somewhere 
in between the idealised series and parallel values (Figure 3.5). More detailed correlations have 
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been developed for some porous structures (Rohsenow ef al., 1998), alternatively this could be 
calculated directly using CFD software if the geometry is explicitly known to validate the porous 
zone properties. 


Hot Cold Hot Cold 
— — od ———— 
(a) Series (Perpendicular) (b) Parallel (Aligned) 
ke = Way eto ker = (1 — p)ks + hke 


Figure 3.5: Bounding values for effective thermal conductivity. 


In some cases, it may be more appropriate to use anisotropic values of thermal conductivity based 
on the geometry (e.g. rod bundle). For simple aligned solids, the thermal conductivity can be es- 
timated using the series/parallel values (i.e. perpendicular to and aligned with the rod bundle re- 
spectively). 


Depending on the porous structure (proportion of solid, open grid or random structure), radiative 
heat transfer is likely to be absorbed by a porous medium. Therefore, in most cases, the porous 
zone can be considered as a solid region that absorbs and emits thermal radiation. 


For thermal radiation within the porous zone, the decomposition (T}—T}) = (T1-T2)(T1 +T2)(T7 + 
T2) may be useful for creating equivalent thermal conductivity models in regions modelled with a 
porous representation where thermal radiation transport occurs, because if the approximation can 
be made that T; ~ T2 = T, then 


ae ee Vas OF 


This allows an equivalent conductivity model to be created that depends on local temperature to the 
third power to represent the contribution of thermal radiation, without modelling thermal radiation 
explicitly. 


3.4.4 LES Aspects 


LES is introduced in Section 3.4.1 and meshing aspects are considered in Section 3.4.2.4. This 
section presents some advice for conducting LES work for CHT analysis with further details on the 
general application of LES modelling provided in Volume 3 (Section 3.2.5). Despite computational 
costs (discussed below) the application of LES to NTH is steadily increasing. 


Unsteady CHT occurs due to the fluctuating heat loads generated by the turbulent flow. Simulations 
of turbulent flows using LES may capture the flow unsteadiness and associated variable heat loads 
on structures. Some industrial examples of LES with CHT include: 


¢ Thermal fatigue and striping in T-junctions (Kuhn et al., 2010, Ferrara and Di Marco, 2017); 
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¢ Round jet impinging on a flat plate (HadZiabdi¢é and Hanjali¢, 2008); 
« Fuel cladding temperatures in wire-wrapped fuel assembly (Doolaard et al., 2015); 
* Thermal mixing in annulus of control-rod guide tube (Bergagio et al., 2020). 


Computational Cost: As noted in Volume 3, the computational cost of LES is less than DNS 
but significantly higher than RANS. This is particularly so for high Re or Gr boundary layer flows, 
and because near-wall turbulence levels are coupled to surface heat transfer, the computational 
cost also increases as Ra increases. For CHT analysis, the thermal mass associated with the solid 
(Section 2.1.1.2) means that the temperature in the solid can take a long time to stabilise. 


Because of the inherent computational cost associated with this method, the simulated physical 
time with such approaches is limited (typically several seconds up to one minute in most engineer- 
ing applications, Koren et al., 2017). Therefore, simulating low-cycle variations and long transients 
using LES is unlikely to be feasible currently. 


Time Step: LES calculations are transient, so a time step is required by the solver. This is de- 
termined by the time step required for the flow field, as discussed in Volume 3 (Section 3.2.5), 
because the timescales associated with the flow are generally significantly smaller than any tem- 
perature variations within the solid. 


Some solvers or coupling approaches allow different time steps to be used in the fluid and solid 
regions. This has the potential to significantly reduce the computational time required, for example 
by solving the solid domain every 50 iterations of the fluid domain. However, since the wall heat 
flux and temperature depends on the size of the coupling time step, a large value might impact the 
accuracy of an LES prediction as the unsteady phenomena may be under-resolved at the fluid-solid 
interface. 


RANS Aspects 


RANS is introduced in Section 3.4.1 and meshing aspects are considered in Section 3.4.2.5. This 
section presents some advice for conducting RANS work for CHT analysis with further details on 
the general application of RANS modelling provided in Volume 3 (Section 3.2.6). 


Steady and Unsteady RANS: One of the key aspects of RANS is the ability to run models as 
either steady-state (RANS) or unsteady (URANS). While turbulent flows are, strictly, always un- 
steady (due to turbulent fluctuations in flow quantities), large-scale unsteady motion can also arise 
from non-turbulence phenomena or be imposed externally by time-varying boundary conditions. 
This large-scale ‘coherent’ unsteadiness can be considered distinct from the unsteady small-scale 
‘incoherent’ unsteadiness typically associated with turbulence. URANS methods exploit this dis- 
tinction by solving a form of the RANS equations which retain the time-derivative term. 


Steady-state calculations account for the majority of engineering analysis, but the turbulence mod- 
els used in RANS generally cannot account for large-scale flow motion. Therefore, a steady RANS 
solution may not give the same results as a long-term time-averaged URANS solution with the 
same steady boundary conditions. 


As noted in Volume 1, for many decades URANS is likely to be the only practical CFD method for 
predicting long duration problems such as long reactor shut-down transients. This is particularly 
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true for CHT analyses as the timescales associated with low-cycle variations and long thermal 
transients are significant and can occur over several hours. 


For flows whose stability is not certain, CHT should be carefully considered. Since the thermal 
properties of the solid have the potential to be substantially different from the fluid, it is worth un- 
derstanding whether any large-scale flow motion has any impact on the solid temperatures and so 
whether a URANS analysis is necessary. In addition, due to the slow thermal response associated 
with solids it is often important to start any URANS analysis from a properly converged steady 
RANS solution where the solid temperature field has also reached steady-state. This may not be 
the appropriate choice in some cases; if modelling a relatively sudden change in plant conditions is 
the purpose of the simulation, then large thermal mass objects that are not in thermal equilibrium 
with either the initial or subsequent flow can significantly affect the response. 


Time Step: If transient data is required from a URANS solution, a time step is normally specified. 
The choice of time step will depend on the timescales present in the flow, and can have a significant 
impact on the results. Estimating a time step based on the flow field is discussed in Volume 3 
(Section 3.2.6). 


As for LES, the timescales associated with the fluid motion are likely to be significantly smaller than 
the timescales associated with the solid. However, this means that the timescales required for the 
solid to reach a statistically steady value will be significantly higher, and so the model may need 
to run for a long time especially if simulating a long thermal transient. Depending on the coupling 
approach, it may be possible to use different time steps in the fluid and solid regions to reduce the 
overall computational time. 


For coupled monolithic solvers, to minimise the time taken to solve each iteration the fluid and solid 
regions should be decomposed for parallel solution so that the computational effort is balanced, 
although CFD codes should be able implement this themselves without user intervention. 


Turbulence and Near-Wall Modelling: RANS turbulence models and near-wall modelling ap- 
proaches are considered in Volume 3 (Sections 2.2.4 and 3.2.6) as they are key aspects for buoy- 
ancy affected flows. Since the presence of a solid region does not fundamentally change the NTH 
phenomena being modelled, the most appropriate turbulence model should be selected based on 
the flow field, although some RANS models are being developed for CHT. For example, Craft et al. 
(2010) describe that predicting the correct details of the temperature response in a wall requires 
an accurate representation of the associated near-wall turbulence processes in the fluid. 


For CHT analysis, the heat transfer at the wall is extremely important, and so a wall resolving 
(low-Re, fine mesh, small y*) near-wall modelling approach is recommended to ensure that the 
near-wall thermal boundary layer is properly resolved. 
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Extracting Heat Transfer Data 


It is often necessary to extract heat transfer data from CFD simulations. This can be either as a 
heat flux, surface temperature or HTC with a local bulk fluid temperature, where: 


q= A(T» _ Tw) 


Heat flux and surface temperature are a direct output of solving the energy equation within a CFD 
analysis and so can be easily extracted and used to pass information in a coupled partitioned ap- 
proach (Section 3.3.2). Therefore, this section considers the methods available to calculate HTCs 
from a CFD solution, which are used to: 


« Understand, interpret and assess CFD predictions by comparison against empirical correla- 
tions and theory. 


¢ Undertake an uncoupled CFD/FEA solution (Section 3.3.1) when the thermal environment is 
changing so slowly that the fluids and solids can be considered to be in a quasi-steady-state, 
or the HTC is effectively independent of the solid temperature (i.e. forced convection). 


As discussed in Section 3.3.2, HTCs are a useful way of exchanging heat transfer information, 
but are dependent on the local bulk temperature which is often non-uniform and not easy to de- 
fine. Therefore, a number of different methods have been developed to extract HTCs. The most 
appropriate approach will depend on the local flow field and objectives of the analysis. 


It is worth noting that for all but the simplest scenarios, there is often uncertainty in the calculated 
value of HTC. This needs to be taken into account in any subsequent FEA models. It is also 
important to ensure that the solid geometry, HTCs and bulk temperatures are consistent between 
models to ensure that the calculated temperatures and heat fluxes are appropriate. 


Reynolds analogy: Under specific conditions (no external pressure gradient, Pr = 1 and uniform 
wall temperature), the momentum flux and heat flux can be shown to be analogous, and the velocity 
and temperature profiles have the same shape. This analogy can be used to estimate an HTC as 
follows: 


Cr _ h 5 es CoTw 

2 PasaGa ~ Us 
Where C; is the skin friction coefficient and U. is the free-stream velocity. Chilton-Coloburn devel- 
oped an empirical correlation based on experimental data for a flat plate and fully developed pipe 
flow (Chilton and Colburn, 1934), which is often called the modified Reynolds analogy’: 


¢ =Srer, n= 


CpTw 


ice Pr < 60 
oo Ff 


This analogy is not appropriate for liquid metal heat transfer due to the low Prandtl number as 
discussed in Volume 5. 


The term StPr2/3 = Nu/(RePr/?) is known as a j-factor; the Chilton-Colburn analogy is also applied to mass transfer, 
where Pr and Nu are replaced by Schmidt and Sherwood numbers respectively (Incropera et al., 2011, Bird et al., 2007, 
Chapter 22). 
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Bulk temperature: For forced convection in a pipe or mixed convection within an enclosed cavity, 
it is often convenient to calculate a ‘bulk’ or ‘mean’ static temperature (T,) based on an appropriate 
region (Section 2.1.2): 


« If the fluid near a surface is well mixed, a single constant value of bulk temperature may be 
appropriate. 


¢ For the flow along a pipe/duct, an offset line along the pipe axis should be used with the bulk 
temperature calculated based on the local average fluid temperature. 


« An offset surface may also be used to determine a local fluid temperature and calculate a 
local HTC. The offset distance should be set so that it is at the edge of the thermal boundary 
layer and is only valid when thermal gradients in the free-stream are small. Since this does 
not take into account any variation in boundary layer thickness, a more general approach 
based on the y* at the edge of the boundary layer is discussed below. 


This is a simplistic approach that needs to be used carefully as it can generate unrealistic or 
negative values of HTC when the local wall temperature equals the bulk temperature. Therefore, it 
is often appropriate to calculate an average value of HTC based on the average bulk temperature. 
This can provide a more representative estimate of the HTC on a surface. 


Diabatic Y-plus Field (DYF): The local fluid temperature for the DYF is defined as the fluid 
temperature at the edge of the thermal boundary layer, where the thermal gradient is zero. This 
can be extracted from a CFD solution by creating a surface in the CFD model that is offset a fixed 
dimensionless distance from the wall i.e. y* of 60 to 250 (Luczynski et al., 2019). A local HTC can 
then be determined using the local fluid and wall temperature values. 


This method is appropriate for forced convection flows, when the flow is clearly defined with no 
separation or stagnation, so that the edge of the thermal boundary layer can be identified. However, 
if there are large thermal gradients in the flow due to mixed convection or local forced convection 
effects, then it can lead to unrealistic or negative HTC values. 


Thermal superposition: If the energy equation near the wall is linear with temperature, which 
is generally true if buoyancy effects and static temperature changes are small, superposition of 
solutions may be used to show that if the wall temperature is changed from 7, to Tz (Verdicchio, 
2001), then HTC can be calculated explicitly as small changes in local wall temperature have a 
negligible effect on the local HTC and local fluid temperature: 


_ U1 — A 


qi = h(Ty — Tp) and q2 = A(T2— Tp), 2. A= ToT 


This means that an estimate of HTC, and hence local fluid temperature, can be obtained from two 
CFD solutions with different wall temperatures. In practice, this can be achieved by calculating a 
baseline or adiabatic surface temperature in a CFD solution, extracting the temperature profile on 
a small zone or surface and applying a small fixed temperature change to it (e.g. 10°C to 50°C) 
and then re-solving the CFD solution with a fixed temperature boundary condition. 


Since this method uses two solutions to derive the HTC, it can be applied when the local fluid 
temperature cannot be clearly defined or there is a complex free-stream flow field. This has been 
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successfully validated for a gas turbine compressor drum (Lewis and Provins, 2004) to extract 
HTCs from a CFD solution and apply them to an FEA model. 


Experimental Methods 


The role of experimental methods in NTH analysis and the value of experimental and modelling 
teams working together closely is introduced in Volume 1 (Section 4.6). This section builds on 
this, discussing the different techniques that might be used to measure temperature, and some of 
the challenges associated with them, while pressure, velocity and flow visualisation measurement 
techniques are discussed in Volume 3 (Natural Convection and Passive Cooling). 


Temperature Measurement 


A wide range of different temperature measurement techniques are available, with different levels 
of accuracy, response times and accessibility requirements. 


Thermocouple: Thermocouples are based on the principle that a temperature-dependent volt- 
age is produced when dissimilar metals are connected together. The five standard base metal 
thermocouples are chromel-constantan (Type E), iron-constantan (Type J), chromel-alumel (Type 
K), nicrosil-nisil (Type N) and copper-constantan (Type T). These are all relatively cheap to man- 
ufacture, but some can become inaccurate over time (Morris and Langari, 2020). This drift in cali- 
bration over time is less of a problem in experiments than in actual reactors, where thermocouples 
are expected to provide readings for decades. In addition, noble metal thermocouples have been 
developed that have high stability and long life even at high temperatures, but they are expensive 
and cannot be used in reducing atmospheres. In most cases, base metal thermocouples are used, 


which have a typical accuracy of +0.5% to + 0.75%, and the most appropriate type will depend 
on the environment, stability, operating life and accuracy requirement. 


It is often impractical to connect a voltage-measuring instrument in close proximity to the point at 
which temperature is being measured, and so extension leads up to several meters long are re- 
quired. This modifies the electrical circuit, and so care must be taken regarding the materials used 
in the extension leads and calibration of the thermocouple. If transient measurements are required 
the tip should be left exposed to maximise the speed of response, although they are delicate de- 
vices that must be treated carefully. Thermocouples can be protected using an insulated sheath 
(often called a probe), although this increases the time taken for the thermocouple to respond (e.g. 
0.15 s time constant for 1mm diameter sheath). 


Thermocouples respond quickly to changes in temperature, and are a simple, cheap and widely 
used method for measuring fluid and solid surface temperature. 


Resistance Thermometer: Resistance thermometers or Resistance Temperature Devices (RTDs) 
are based on the principle that the resistance of a metal varies with temperature. Platinum has a 
good linear resistance/temperature characteristic and it also has good chemical inertness (Mor- 
ris and Langari, 2020). In the case of non-corrosive and non-conducting environments, RTDs are 
used without protection. In all other applications, they are protected inside a sheath, which reduces 
the speed of response to changes in temperature. 
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RTDs are commonly used up to 650°C with a measurement inaccuracy of +0.5%. They are 
generally more stable than thermocouples and can measure small temperature differences. 


Distributed Temperature Sensor: Distributed Temperature Sensors (DTSs) use fibre optic ca- 
bles to measure the temperature profile along the sensor cable because changes in temperature 
along the length of the cable alter the amount of light that is transmitted. Analysis of the backscat- 
tered radiation then enables temperature versus distance to be determined (Morris and Langari, 
2020). High resolution DTSs have a high measurement accuracy, and have been used to mea- 
sure temperature distributions near the core of a research reactor and in corrosive liquid-sodium 
environments (Gerardi et al., 2017). 


Thermography: Thermography, or thermal imaging, involves scanning an IR radiation detector 
across an object. This is a non-invasive technique used to provide temperature measurements at 
a single point or on a surface area (Kim et a/., 2016). The IR detector uses the principle that all 
objects emit electromagnetic radiation as a function of temperature. The temperature at the point 
that the instrument is focused on is inferred from a measurement of the incoming infrared radiation. 


Measurement resolution can be high, with temperature differences as small as 0.1°C being de- 
tectable (Morris and Langari, 2020). However, the thermal radiation from a body is very sensitive 
to the composition and surface condition of the body (i.e. its emissivity). Therefore, all IR detectors 
have to be carefully calibrated for each particular body or by applying patches of known emissivity. 
In addition, IR thermography requires line of sight with the emitting body and the results can be 
affected by absorption and scattering between the emitting body and the detector (e.g. dust and 
water droplets). 


IR thermography has high sensitivity and fast response time, and can therefore be exploited to 
effectively measure convective heat fluxes with steady or transient techniques, or to perform de- 
tailed thermal surface flow visualisation (Astarita and Carlomagno, 2013). IR thermography can 
also be used in situations where contact measurements would be hazardous or impossible (e.g. 
at temperatures exceeding 1300 °C, Kim ef a/., 2016). However, the method is relatively expensive 
and can only be used to measure the temperature on a solid surface with gas present between the 
surface and infrared detector. 


Temperature Sensitive Paint: Temperature Sensitive Paint (TSP) is a method for measuring 
surface temperature, usually in aerodynamic settings (Liu, 2011). TSP is a paint-like coating that 
fluoresces under a specific illumination wavelength in differing intensities depending on the sur- 
face temperature due to thermal quenching that reduces the luminescent intensity. TSP is able 
to provide non-contact, high resolution, quantitative mapping of surface temperature on complex 
surfaces at low cost. A generic TSP measurement system is composed of paint, illumination light, 
photodetector and data acquisition/processing unit, with each pixel of the imaging camera acting 
like an individual thermocouple. 


Time-resolved TSP applications involve pulsed excitation, and so depends on the ability to effec- 
tively detect the emission from TSP in short exposures. One can thus determine surface temper- 
ature variation as a function of time. In this case, the imaging devices must be synchronized to 
the excitation with multi-channel digital delay/pulse generators providing that synchronization. This 
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technique is used in both aerodynamic and turbomachinery applications, and can be used to iden- 
tify the transition from laminar to turbulent flow due to the abrupt change from low to high heat 
transfer. 


Since TSP uses imaging cameras, the test fluid and test section walls must be optically transparent, 
which makes it difficult to use at reactor pressures and temperatures. 


Measurement Challenges 


Use of test rigs to measure temperature and CHT effects at a system and component level can 
be challenging, as the material properties of the structure around the test section is a fundamental 
part of the test requirement. This limits the accessibility and visibility of the flow within the test 
section, making detailed concurrent flow and thermal measurements difficult. 


If the transient response of the solid is important (high frequency fluctuations or long fault scenario), 
then the scaling of the rig is further complicated by the need to match the solid thermal mass as 
well as the flow phenomena. Other factors include understanding and minimising the impact of 
the external ambient conditions (e.g. heat loss), measurement accuracy and ensuring that the 
instrumentation does not impact the test results. 


Thermocouples are often used within test rigs as they are simple to use and provide point transient 
temperature measurements of the fluid and solid, but care needs to be taken when installing and 
calibrating the instrumentation to ensure that: 


* The thermocouple is in good contact with the surface, otherwise it can be influenced or offset 
by the fluid temperature. Surface measurements can be particularly difficult as they can be 
affected by the local flow and small gaps can form over time. 


* The thermocouple does not influence the flow/temperature as they are about 1 mm diameter 
and attached to a cable. 


* The responsiveness (time constant) of the thermocouple is appropriate for the timescales of 
temperature fluctuations being measured. 


* The thermocouples are calibrated and the uncertainty is quantified, as this can be significant 
compared to the temperature changes being measured (approx. 2°C). 


* The location of the thermocouples is optimised to identify and measure the expected features 
of interest in the test. 


As discussed in Volume 1, the most important aspect of a test rig designed for model validation is to 
ensure that the experimentalists and model developers work closely together throughout the design 
process to ensure that the rig meets the model validation requirements and that uncertainties are 
quantified and documented. 
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To conclude the present volume, this section considers how the aspects presented on the preced- 
ing pages may change in the future. 


It is expected that the ongoing exponential increase in computational power and CFD software 
development will enable more accurate CFD predictions of convection, thermal radiation and con- 
jugate heat transfer. As a result, there will be increased use of CFD methods to reduce cost, 
increase performance and improve the safety of advanced reactor designs. This is likely to include: 


Larger (amount of plant modelled) and longer (length of transient) RANS/URANS simula- 
tions. 


More widespread use of LES and hybrid methods, including to generate high-fidelity results 
and benchmark solutions for RANS/URANS model development and validation (Volume 3, 
Section 4). 


Automatic geometry and mesh generation for complex geometries to reduce the cost and 
timescales associated with CFD methods. 


Improved fluid-solid coupling to make it easier to include solid geometry and solve long ther- 
mal transients. 


However, NPPs are complex, integrated systems involving a wide range of different physics, timescales 


and fault scenarios that need to be designed and assessed at a system level. Although CFD meth- 
ods offer significant benefit in terms of accuracy and confidence, they are unlikely to be able to 
routinely simulate all of the details of a whole primary circuit in the near future. Therefore, the in- 
tegrated digital design process and coarse grid CFD should be points of focus for future research 
intended to support the design of advanced reactors. 


Integrated Digital Design Process 


As discussed in Section 3.1, multiscale NTH analysis and multiphysics models that include thermal 
hydraulics are being increasingly developed and used as part of an integrated framework for reactor 
design. This has the potential to simulate the whole primary circuit with CFD resolution in local 
areas, and include coupled neutronics, structure and chemistry (CRUD, corrosion) calculations. 
Integrated frameworks that are currently being developed include: 


MOOSE: Multiphysics Object-Oriented Simulation Environment (MOOSE)' is a multiphysics frame- 
work that has been developed by Idaho National Laboratory (INL) to provide a high-level in- 
terface for computational analysis. It provides a fully-coupled, fully-implicit solver that allows 


1 moose.inl.gov/SitePages/Home.aspx 
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independent codes to be coupled and exchange information, and has the capability to solve 
problems involving solid mechanics, phase field modelling, heat conduction and fluid flow. 

It includes a number of thermal hydraulic analysis software, such as RELAP-7, SAM, TRACE 
and Nek5000, which are described in Volume 1 to enable multiscale analyses. In addition, the 
United States Nuclear Regulatory Commission (US NRC) has proposed a code suite for non- 
LWR systems safety analysis (US NRC, 2020), called the Comprehensive Reactor Analysis 
Bundle (CRAB), which is based around MOOSE. This documents the US NRC’s vision and 
strategy for safety analysis and licensing of Advanced Nuclear Technologies (ANTs). 


SCALE SERPENT FLUENT 
Cross-sections Cross-sections CFD 
PARCS MAMMOTH PRONGHORN Nek5000 
Neutronics Neutronics Core T/H CFD 
TRACE MOOSE 
System and Core T/H Coupling, Tensor Mechanics 
FAST |... =! BISON SAM 
Fuel Performance Fuel Performance System and Core T/H 


Planned Completed 


NRC Code Int'l Code DOE Code Coupling Coupling Input/BC Data 
eS 


oe > 


Figure 4.1: US NRC Comprehensive Reactor Analysis Bundle (CRAB) (US NRC, 2020). 
This image is reproduced courtesy of US NRC, and is not covered by the creative com- 
mons license defined in the legal statement for the present document. 


NURESAFE: NUclear REactor SAFEty simulation platform (NURESAFE) was a European Com- 

mission (EC) funded collaborative research project that was completed in 2015 and in- 
volved 22 organisations (Chanaron, 2016). This built on previous EC Framework Programs 
(NURESIM and NURISP) and was intended to enhance the prediction capability for the de- 
sign and safety assessment of LWR NPPs through the dynamic 3D coupling of multiphysics 
codes into a common simulation scheme. 
The NURESIM platform, at the end of the NURESAFE project, includes neutronics, core 
simulators and fuel thermo-mechanics, as well as system (CATHARE, ATHLET), subchan- 
nel (CTF, SUBCHANFLOW) and CFD (NEPTUNE_CFD, TRIO_U) codes. This is based 
on the open source simulation platform SALOME? with uncertainty quantification analysis 
(URANIE) and includes a generic coupling interface called ICOCO (Interface for Code Cou- 
pling). Although this project was specifically focused on simulating normal operation and 
design basis accidents for LWRs, many of the codes within the platform are being developed 
and validated for ANTs. 


ONCORE: The International Atomic Energy Agency (IAEA) Open-source Nuclear Codes for Reac- 
tor Analysis (ONCORE)’ initiative is an international collaboration framework for the develop- 


2 www.salome-platform.org 
3 nucleus.iaea.org/sites/oncore/SitePages/Home.aspx 
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ment and application of open source multiphysics simulation tools. This includes neutronics, 
structural mechanics, fuel behaviour and thermal hydraulics analysis tools. ONCORE offers 
access to knowledge and tools that otherwise may not easily be available, and promotes 
collaboration and sharing of resources. 


The objective is to develop a consistent open source platform with verification through stan- 
dard problems and code-to-code comparison and documentation. It includes application 
frameworks for integrating multiple codes and a number of the codes are based on the Open- 
FOAM software package. 


Further integration and cost benefit could be achieved in the future through an Integrated Nuclear 
Digital Environment (INDE) that covers the whole life-cycle of the reactor from prototype design 
through operations and decommissioning to storage and waste disposal (Patterson et a/., 2016). 
This would enable the model fidelity to increase as the reactor design develops and track the 
evolution of the geometry. In addition, multiscale and multiphysics models could be linked to real- 
world data from validation of prototypes, in-service monitoring and plant inspections to produce a 
digital twin of the reactor. 


For example, the Nuclear Virtual Engineering Capability (NVEC)* project funded by the UK Gov- 
ernment Department for Business, Energy and Industrial Strategy (BEIS) Nuclear Innovation Pro- 
gramme (NIP) (Volume 1, Appendix A) is using real project applications to demonstrate improved 
efficiency and position the UK as a leader in digital twin and virtual engineering. 


Coarse Grid CFD 


It is currently not feasible to use CFD to perform core-level design calculations due to its pro- 
hibitively high computational cost. Therefore, various solutions have been proposed, among which 
the low resolution approaches show great potential, because they not only retain the flexibility and 
3D advantages of conventional CFD, but also balance the simulation accuracy and computing cost. 
These low resolution approaches can be roughly divided into two categories: 


Porous media method: This uses porous media with uniform or non-uniform porosity to model 
the fluid and solid domain as a whole without representing the detailed geometry. The in- 
bundle friction, spacer induced mixing, and fluid-solid heat transfer are all modelled using 
empirical correlations. Under certain situations, this method can be seen as a flexible sub- 
channel approach implemented on the CFD framework (Gerschenfeld et a/., 2019, Yoon 
et al., 2018). 


Reduced resolution CFD: The fuel assembly geometry is resolved using a coarse body-fitted 
grid. To overcome the shortcomings of the coarse grid in capturing the near-wall effects, 
appropriate modelling strategies, empirical near-wall closure methods or other special treat- 
ments are usually required. The method proposed by Glass ef al. (2011), referred to as 
Coarse-Grid CFD (CG-CFD) uses as few as 5 cells per sub-channel cross-section. Instead 
of solving the Navier-stokes governing equation, they solve an Euler equation with a vol- 
umetric force term implemented to mimic the effects of friction and turbulence. However, 
the derivation of this volume force relies on a large number of pre-performed resolved CFD 


4 n-vec.co.uk 
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simulations on the same geometry (or a representative geometry) and a subsequent pa- 
rameterisation process of the flow conditions. Hanna et al. (2020) introduced an advanced 
machine learning approach to predict CG-CFD simulation errors, though work is at an initial 
stage and a long way from industry application. 


An alternative approach to reduced resolution CFD is to use a novel CFD-based subchannel frame- 
work (Liu et al/., 2019). This combines the advantages of modern CFD and traditional subchannel 
codes, and has been developed as part of the BEIS thermal hydraulics NIP (Volume 1, Appendix 
A)°. This method includes the capability to embed resolved CFD models when local refinement is 
required, but is implemented in a single code to simplify the spatial and temporal data transfer. 


Since this method uses a very coarse mesh, compared to traditional CFD approaches, large re- 
actor components or even the whole reactor core can potentially be modelled. The wall modelling 
closure in Subchannel CFD (or ‘SubChCFD’) is achieved using well validated industry standard 
correlations, similar to those used in subchannel codes, to improve the reliability of the prediction 
without using a fine near-wall mesh. SuobChCFD uses a two-level mesh system (Figure 4.2): 


1. The filtering mesh, which is equivalent to the mesh of a traditional subchannel code, enabling 
the use of existing engineering correlations to account for the integral wall shear and heat 
transfer effects. 


2. The computing mesh, on which a typical CFD solver can be used to resolve the inviscid flow 
with corrections for diffusion using a simple mixing length turbulence model, which can be 
calibrated for specific applications. 


Figure 4.2: Mesh system in SubChCFD. 


SubChCFD has been implemented in the Code_Saturne open source CFD code (but could be 
implemented as an add-on module in any CFD tool) and validated against a number of fuel bun- 
dle test cases. It is able to predict both the velocity and temperature fields for both standard and 
blocked fuel assemblies (Liu, 2019). The ability to include embedded resolved CFD models and 


5 www.innovationfornuclear.co.uk/nuclearthermalhydraulics.html 
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porous media has recently been implemented and validated to locally improve the simulation reso- 
lution where the flow and heat transfer exhibit complex features (Liu, 2021). This method is showing 
significant promise, and in the future could provide a more efficient and accurate way to model the 
reactor core and large reactor systems within a unified CFD framework. 
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6 Nomenclature 


Latin Symbols 
A Area, m? 
At Atwood number (At = (p1 — p2)/(e1 + p2)) 
Bi Biot number (Bi = hL/kg) 
Cp, Cy Specific heat at constant pressure or volume, J kg~! K~* 
dorD_ Diameter (D, = 4A;;/Pcs for hydraulic diameter), m 
f Darcy-Weisbach friction factor 
Fo Fourier number (Fo = a,t/L?) 
Gr Grashof number (Gr = gL3Ap/v?p = gL?BAT/v?, using the Boussinesq approxi- 
mation Ap/p ~ —BAT, where AT is often taken as Ty, — Ts,00) 
g Acceleration due to gravity, ms~? 
h_ Specific enthalpy, J kg~+, Heat Transfer Coefficient (HTC), Wm~? K—! or height, m 
I Radiative intensity, W m~? sr~! or Wm~? sr~+ wm? for a spectral density, where sr 
(steradian) is solid angle 
J Radiosity, W m~? 
k Thermal conductivity, Wm-! K-? 
L Length or wall thickness, m 
M_ Molar mass of a species, kg kmol~? 
Ma _ Mach number (Ma = U/a, where a is the speed of sound) 
n_ Refractive index 
Nu Nusselt Number (Nu = AL/k¢) 
p Perimeter, m 
P Pressure (P; = static pressure, Py = total pressure), Nm~? or Pa 
Pe  Péclet number (Pe = RePr) 
Pr Prandtl number (Pr = cpu/ ke) 
q_ Heat flux (rate of heat transfer per unit area, g = Q/A), Wm~? 
Q_ Rate of heat transfer, W 
r Radius, m 
R_ Gas constant (for a particular gas, R = R/M), Jkg~!K~! 
R Universal gas constant, 8314.5 J kmol~+ K~+ 
Rip Thermal resistance, K W~! 
Ra Rayleigh number (Ra = GrPr) 
Re Reynolds number (Re = pUL/,, or for an internal flow Re = WD,,/A.<h) 
Ri Richardson number (Ri = Gr/Re’) 
Sr Strouhal number (Sr = fL/U, where f is frequency) 
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Nomenclature ¢, 


«x j<< 


Volume 2 


Stanton number (St = Nu/RePr) 

Time, s 

Temperature (7, = static temperature, 77 = total temperature), K 
Wall friction velocity (u, = \/Tw/p), ms~+ 

Velocity, ms~? or thermal transmittance, Wm? K~! 
Specific volume, m* kg~? 

Volume, m? 

Mass flow rate, kg s~+ 
Wall distance, m 


Non-dimensional wall distance (y* = yu, /v) 


Greek Symbols 


Qn SS FE YM KF Hm R WR 


4 


oa 


Thermal diffusivity (a = k/pcp), m?s~? 


Volumetric thermal expansion coefficient (8 = —(1/)(0p/0T)), K~* 
Ratio of specific heats (y = cp/c,) 

Emissivity or surface roughness height, m 

Absorption coefficient, m~? 

Wavelength, m 


Viscosity, kg m~1s~? 


Kinematic viscosity and momentum diffusivity (v = w/p), m?s~? 
Density, kg m~? 

Stefan Boltzmann constant, 5.67 x 10-8 W m-? K~* 

Shear stress, Nm~? 


Porosity or void fraction 


Subscripts and Modifications 


b 
cs 
f 


i 


wn 


8 


~ 


Bulk (mass-averaged) quantity 
Cross-sectional quantity 

Quantity relating to a fluid 

Quantity relating to a particular species 
Total (stagnation) quantity 

Turbulent quantity 

Static quantity or quantity relating to a solid 
Quantity relating to a wall or surface 
Quantity far from a wall or in free-stream 
Average quantity 

Molar quantity 

Varying quantity 
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7 Abbreviations 


AGR 
ANT 
AOA 
BEIS 
BWR 
CAD 
CASL 
CCFL 
CFD 
CG-CFD 
CHT 
CIPS 
CRAB 
CRUD 
CSNI 
DNS 
DO 
DTS 
DYF 
EC 
ECCS 
FAC 
FCH 
FEA 
FIV 
FSI 
HTC 
IAEA 


Advanced Gas-cooled Reactor 

Advanced Nuclear Technology 

Axial Off-set Anomaly 

Department for Business, Energy and Industrial Strategy 
Boiling Water Reactor 

Computer Aided Design 

Consortium for Advanced Simulation of Light Water Reactors 
Counter-Current Flow Limitation 

Computational Fluid Dynamics 

Coarse-Grid CFD 

Conjugate Heat Transfer 

CRUD Induced Power Shift 

Comprehensive Reactor Analysis Bundle 

Chalk River Unidentified Deposits 

Committee on the Safety of Nuclear Installations 
Direct Numerical Simulation 

Discrete Ordinates 

Distributed Temperature Sensor 

Diabatic Y-plus Field 

European Commission 

Emergency Core Cooling System 
Flow-Accelerated Corrosion 

First Cell Height 

Finite Element Analysis 

Flow Induced Vibration 

Fluid-Structure Interaction 

Heat Transfer Coefficient 

International Atomic Energy Agency 

Integral Effect Test 

Integrated Nuclear Digital Environment 

Idaho National Laboratory 

Infrared 

International Organization for Standardization 
Lead-Bismuth Eutectic 

Large Eddy Simulation 

Liquid Metal-cooled Fast Reactor 
Loss-Of-Coolant Accident 

Light Water Reactor 

Monte Carlo 

Multiphysics Object-Oriented Simulation Environment 
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Abbreviations Oo 


MSR 
NDT 

NIP 

NPP 
NTH 
NURESAFE 
NVEC 
ONCORE 
PCT 
PHRS 
PTS 
PWR 
RANS 
RMI 
RMS 
ROM 
RPV 
RTD 
RTE 
S28 
SFR 
SGS 
SSC 
TSP 
URANS 
US NRC 
V&V 
VERA 
WALT 
WMLES 


Volume 2 


Molten Salt Reactor 

Non-Destructive Testing 

Nuclear Innovation Programme 

Nuclear Power Plant 

Nuclear Thermal Hydraulics 

NUclear REactor SAFEty simulation platform 
Nuclear Virtual Engineering Capability 
Open-source Nuclear Codes for Reactor Analysis 
Peak Cladding Temperature 

Passive Decay Heat Removal System 
Pressurised Thermal Shock 

Pressurised Water Reactor 
Reynolds-Averaged Navier-Stokes 
Reflective Metallic Insulation 

Root Mean Square 

Reduced Order Model 

Reactor Pressure Vessel 

Resistance Temperature Device 

Radiative Transfer Equation 
Surface-to-Surface 

Sodium-cooled Fast Reactor 
Sub-Grid-Scale 

Structure, System and Component 
Temperature Sensitive Paint 

Unsteady Reynolds-Averaged Navier-Stokes 
United States Nuclear Regulatory Commission 
Verification and Validation 

Virtual Environment for Reactor Applications 
Westinghouse Advanced Loop Tester 

Wall Modeled Large Eddy Simulation 
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